25 #if defined (HAVE_CONFIG_H)
93 operator R () {
return scl *
std::pow (sum, 1/p); }
105 template <
typename U>
121 operator R () {
return scl *
std::pow (sum, -1/p); }
125 template <
typename R>
154 operator R () {
return scl * std::sqrt (sum); }
158 template <
typename R>
164 template <
typename U>
169 operator R () {
return sum; }
173 template <
typename R>
179 template <
typename U>
187 operator R () {
return max; }
191 template <
typename R>
197 template <
typename U>
205 operator R () {
return min; }
209 template <
typename R>
215 template <
typename U>
218 if (val != static_cast<U> (0)) ++
num;
220 operator R () {
return num; }
225 template <
typename T,
typename R,
typename ACC>
235 template <
typename T,
typename R,
typename ACC>
243 accj.accum (m(
i, j));
245 res.
xelem (j) = accj;
249 template <
typename T,
typename R,
typename ACC>
253 std::vector<ACC> acci (m.
rows (), acc);
257 acci[
i].accum (
m(
i, j));
265 template <
typename T,
typename R,
typename ACC>
273 accj.accum (m.
data (
k));
275 res.
xelem (j) = accj;
279 template <
typename T,
typename R,
typename ACC>
283 std::vector<ACC> acci (m.
rows (), acc);
295 #define DEFINE_DISPATCHER(FUNC_NAME, ARG_TYPE, RES_TYPE) \
296 template <typename T, typename R> \
297 RES_TYPE FUNC_NAME (const ARG_TYPE& v, R p) \
301 FUNC_NAME (v, res, norm_accumulator_2<R> ()); \
303 FUNC_NAME (v, res, norm_accumulator_1<R> ()); \
304 else if (lo_ieee_isinf (p)) \
307 FUNC_NAME (v, res, norm_accumulator_inf<R> ()); \
309 FUNC_NAME (v, res, norm_accumulator_minf<R> ()); \
312 FUNC_NAME (v, res, norm_accumulator_0<R> ()); \
314 FUNC_NAME (v, res, norm_accumulator_p<R> (p)); \
316 FUNC_NAME (v, res, norm_accumulator_mp<R> (p)); \
329 template <typename ColVectorT, typename R>
338 R
fi =
i *
static_cast<R
> (M_PI) / nsamp;
339 R lambda1 = cos (fi);
343 lambda1 /= lmnr; mu1 /= lmnr;
357 template <
typename ColVectorT,
typename R>
361 std::complex<R>& lambda, std::complex<R>& mu)
363 typedef std::complex<R> CR;
366 CR lamcu = lambda /
std::abs (lambda);
371 R
fi =
i *
static_cast<R
> (M_PI) / nsamp;
372 R lambda1 = cos (fi);
376 lambda1 /= lmnr; mu1 /= lmnr;
377 R nrm1 =
vector_norm (lambda1 * lamcu * y + mu1 * col, p);
380 lambda = lambda1 * lamcu;
390 R
fi =
i *
static_cast<R
> (M_PI) / nsamp;
391 lamcu = CR (cos (fi), sin (fi));
392 R nrm1 =
vector_norm (lama * lamcu * y + mu * col, p);
395 lambda = lama * lamcu;
402 template <
typename T,
typename R>
412 template <
typename VectorT,
typename R>
415 VectorT res (x.dims ());
422 template <
typename MatrixT,
typename VectorT,
typename R>
423 R
higham (
const MatrixT&
m, R
p, R tol,
int maxiter,
426 x.resize (m.columns (), 1);
428 VectorT
y(m.rows (), 1, 0), z(m.rows (), 1);
429 typedef typename VectorT::element_type RR;
435 VectorT col (m.column (
k));
441 y = lambda *
y + mu * col;
450 while (iter < maxiter)
461 || (gamma - gamma1) <= tol*gamma))
482 template <
typename MatrixT,
typename VectorT,
typename R>
499 const R sqrteps = std::sqrt (std::numeric_limits<R>::epsilon ());
500 res =
higham (m, p, sqrteps, max_norm_iter, x);
509 template <
typename MatrixT,
typename VectorT,
typename R>
520 const R sqrteps = std::sqrt (std::numeric_limits<R>::epsilon ());
521 res =
higham (m, p, sqrteps, max_norm_iter, x);
531 #define DEFINE_XNORM_FUNCS(PREFIX, RTYPE) \
532 OCTAVE_API RTYPE xnorm (const PREFIX##ColumnVector& x, RTYPE p) \
534 return vector_norm (x, p); \
536 OCTAVE_API RTYPE xnorm (const PREFIX##RowVector& x, RTYPE p) \
538 return vector_norm (x, p); \
540 OCTAVE_API RTYPE xnorm (const PREFIX##Matrix& x, RTYPE p) \
542 return svd_matrix_norm (x, p, PREFIX##Matrix ()); \
544 OCTAVE_API RTYPE xfrobnorm (const PREFIX##Matrix& x) \
546 return vector_norm (x, static_cast<RTYPE> (2)); \
555 template <typename T, typename R>
565 #define DEFINE_XNORM_SPARSE_FUNCS(PREFIX, RTYPE) \
566 OCTAVE_API RTYPE xnorm (const Sparse##PREFIX##Matrix& x, RTYPE p) \
568 return matrix_norm (x, p, PREFIX##Matrix ()); \
570 OCTAVE_API RTYPE xfrobnorm (const Sparse##PREFIX##Matrix& x) \
573 array_norm_2 (x.data (), x.nnz (), res); \
580 #define DEFINE_COLROW_NORM_FUNCS(PREFIX, RPREFIX, RTYPE) \
581 extern OCTAVE_API RPREFIX##RowVector \
582 xcolnorms (const PREFIX##Matrix& m, RTYPE p) \
584 return column_norms (m, p); \
586 extern OCTAVE_API RPREFIX##ColumnVector \
587 xrownorms (const PREFIX##Matrix& m, RTYPE p) \
589 return row_norms (m, p); \
OCTAVE_API ColumnVector xrownorms(const Matrix &m, double p)
#define DEFINE_XNORM_FUNCS(PREFIX, RTYPE)
Octave interface to the compression and uncompression libraries.
octave_idx_type rows(void) const
#define DEFINE_DISPATCHER(FUNC_NAME, ARG_TYPE, RES_TYPE)
octave_idx_type numel(void) const
Number of elements in the array.
identity matrix If supplied two scalar respectively For allows like xample val
VectorT dual_p(const VectorT &x, R p, R q)
OCTAVE_NORETURN liboctave_error_handler current_liboctave_error_handler
octave_idx_type * cidx(void)
octave_idx_type columns(void) const
OCTAVE_EXPORT octave_value_list return the number of command line arguments passed to Octave If called with the optional argument the function t
Template for N-dimensional array classes with like-type math operators.
octave_idx_type rows(void) const
static void higham_subp(const ColVectorT &y, const ColVectorT &col, octave_idx_type nsamp, R p, R &lambda, R &mu)
#define DEFINE_COLROW_NORM_FUNCS(PREFIX, RPREFIX, RTYPE)
nd deftypefn *octave_map m
OCTAVE_API RowVector xcolnorms(const Matrix &m, double p)
DM_T singular_values(void) const
octave_int< T > pow(const octave_int< T > &a, const octave_int< T > &b)
R matrix_norm(const MatrixT &m, R p, VectorT)
static const char * p_less1_gripe
void accum(std::complex< R > val)
T & xelem(octave_idx_type n)
N Dimensional Array with copy-on-write semantics.
charNDArray max(char d, const charNDArray &m)
octave_idx_type * ridx(void)
#define DEFINE_XNORM_SPARSE_FUNCS(PREFIX, RTYPE)
=val(i)}if ode{val(i)}occurs in table i
norm_accumulator_mp(R pp)
OCTAVE_EXPORT octave_value_list return the value of the option it must match the dimension of the state and the relative tolerance must also be a vector of the same length tem it must match the dimension of the state and the absolute tolerance must also be a vector of the same length The local error test applied at each integration step is xample roup abs(local error in x(i))<
void array_norm_2(const T *v, octave_idx_type n, R &res)
void row_norms(const MArray< T > &m, MArray< R > &res, ACC acc)
the element is set to zero In other the statement xample y
OCTAVE_EXPORT octave_value_list or N dimensional array whose elements are all equal to the IEEE symbol NaN(Not a Number).NaN is the result of operations which do not produce a well defined 0 result.Common operations which produce a NaN are arithmetic with infinity ex($\infty-\infty $)
void vector_norm(const Array< T > &v, R &res, ACC acc)
R svd_matrix_norm(const MatrixT &m, R p, VectorT)
void column_norms(const MArray< T > &m, MArray< R > &res, ACC acc)
std::complex< float > FloatComplex
R higham(const MatrixT &m, R p, R tol, int maxiter, VectorT &x)
std::complex< double > Complex
Vector representing the dimensions (size) of an Array.
octave_idx_type columns(void) const
F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T const F77_REAL const F77_REAL F77_REAL &F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE &F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T F77_DBLE &F77_RET_T F77_REAL &F77_RET_T F77_REAL &F77_RET_T F77_DBLE &F77_RET_T const F77_DBLE F77_DBLE &F77_RET_T const F77_REAL F77_REAL &F77_RET_T F77_REAL F77_REAL &F77_RET_T F77_DBLE F77_DBLE &F77_RET_T const F77_DBLE * x
charNDArray min(char d, const charNDArray &m)