NcmNNLS

NcmNNLS — Non-negative linear least-squares

Functions

Properties

guint ncols Read / Write / Construct Only
guint nrows Read / Write / Construct Only
double reltol Read / Write / Construct
NcmNNLSUMethod umethod Read / Write / Construct

Types and Values

Object Hierarchy

    GEnum
    ╰── NcmNNLSUMethod
    GObject
    ╰── NcmNNLS

Description

NcmNNLS is a class that implements the non-negative linear least-squares algorithm.

Functions

ncm_nnls_new ()

NcmNNLS *
ncm_nnls_new (guint nrows,
              guint ncols);

Creates a new NcmNNLS object.

Parameters

nrows

number of rows

 

ncols

number of columns

 

Returns

a new NcmNNLS.


ncm_nnls_ref ()

NcmNNLS *
ncm_nnls_ref (NcmNNLS *nnls);

Increase the reference of nnls by one.

Parameters

nnls

a NcmNNLS

 

Returns

nnls .

[transfer full]


ncm_nnls_free ()

void
ncm_nnls_free (NcmNNLS *nnls);

Decrease the reference count of nnls by one.

Parameters

nnls

a NcmNNLS

 

ncm_nnls_clear ()

void
ncm_nnls_clear (NcmNNLS **nnls);

Decrease the reference count of nnls by one, and sets the pointer *nnls to NULL.

Parameters

nnls

a NcmNNLS

 

ncm_nnls_set_umethod ()

void
ncm_nnls_set_umethod (NcmNNLS *nnls,
                      NcmNNLSUMethod umethod);

Sets which unconstrained least-squares method to use.

Parameters

nnls

a NcmNNLS

 

umethod

a NcmNNLSUMethod

 

ncm_nnls_get_umethod ()

NcmNNLSUMethod
ncm_nnls_get_umethod (NcmNNLS *nnls);

Gets the unconstrained least-squares method being used.

Parameters

nnls

a NcmNNLS

 

Returns

a NcmNNLSUMethod.


ncm_nnls_set_reltol ()

void
ncm_nnls_set_reltol (NcmNNLS *nnls,
                     const gdouble reltol);

Sets relative tolerance to reltol .

Parameters

nnls

a NcmNNLS

 

reltol

a double

 

ncm_nnls_get_reltol ()

gdouble
ncm_nnls_get_reltol (NcmNNLS *nnls);

Gets the relative tolerance being used.

Parameters

nnls

a NcmNNLS

 

Returns

the current relative tolerance.


ncm_nnls_get_nrows ()

guint
ncm_nnls_get_nrows (NcmNNLS *nnls);

Parameters

nnls

a NcmNNLS

 

Returns

number of rows.


ncm_nnls_get_ncols ()

guint
ncm_nnls_get_ncols (NcmNNLS *nnls);

Parameters

nnls

a NcmNNLS

 

Returns

number of rows.


ncm_nnls_solve ()

gdouble
ncm_nnls_solve (NcmNNLS *nnls,
                NcmMatrix *A,
                NcmVector *x,
                NcmVector *f);

Solves the system $A\vec{x} = \vec{f}$ for $\vec{x}$ imposing the non negativity constraint on $\vec{x}$, i.e., $\vec{x} > 0$.

Parameters

nnls

a NcmNNLS

 

A

a NcmMatrix $A$

 

x

a NcmVector $\vec{x}$

 

f

a NcmVector $\vec{f}$

 

Returns

the Euclidean norm of the residuals.


ncm_nnls_solve_LH ()

gdouble
ncm_nnls_solve_LH (NcmNNLS *nnls,
                   NcmMatrix *A,
                   NcmVector *x,
                   NcmVector *f);

Solves the system $A\vec{x} = \vec{f}$ for $\vec{x}$ imposing the non negativity constraint on $\vec{x}$, i.e., $\vec{x} > 0$. This method solves the system using the original code by Charles L. Lawson and Richard J. Hanson translated to C using f2c.

Parameters

nnls

a NcmNNLS

 

A

a NcmMatrix $A$

 

x

a NcmVector $\vec{x}$

 

f

a NcmVector $\vec{f}$

 

Returns

the Euclidean norm of the residuals.


ncm_nnls_solve_lowrankqp ()

gdouble
ncm_nnls_solve_lowrankqp (NcmNNLS *nnls,
                          NcmMatrix *A,
                          NcmVector *x,
                          NcmVector *f);

Solves the system $A\vec{x} = \vec{f}$ for $\vec{x}$ imposing the non negativity constraint on $\vec{x}$, i.e., $\vec{x} > 0$. This method solves the system using the LowRankQP quadratic programming code.

Parameters

nnls

a NcmNNLS

 

A

a NcmMatrix $A$

 

x

a NcmVector $\vec{x}$

 

f

a NcmVector $\vec{f}$

 

Returns

the Euclidean norm of the residuals.


ncm_nnls_solve_splx ()

gdouble
ncm_nnls_solve_splx (NcmNNLS *nnls,
                     NcmMatrix *A,
                     NcmVector *x,
                     NcmVector *f);

Solves the system $A\vec{x} = \vec{f}$ for $\vec{x}$ imposing the non negativity constraint on $\vec{x}$, i.e., $\vec{x} > 0$. This method solves the system using function libqp_splx_solver from libqp.

Parameters

nnls

a NcmNNLS

 

A

a NcmMatrix $A$

 

x

a NcmVector $\vec{x}$

 

f

a NcmVector $\vec{f}$

 

Returns

the Euclidean norm of the residuals.


ncm_nnls_solve_gsmo ()

gdouble
ncm_nnls_solve_gsmo (NcmNNLS *nnls,
                     NcmMatrix *A,
                     NcmVector *x,
                     NcmVector *f);

ncm_nnls_get_residuals ()

NcmVector *
ncm_nnls_get_residuals (NcmNNLS *nnls);

Gets the solution residuals, this method return the last residuals computed during ncm_nnls_solve(). If ncm_nnls_solve() was not called the return is undefined.

Parameters

nnls

a NcmNNLS

 

Returns

residuals vector.

[transfer none]

Types and Values

NCM_TYPE_NNLS

#define NCM_TYPE_NNLS (ncm_nnls_get_type ())

enum NcmNNLSUMethod

Method used to solve the intermediate unconstrained least-squares.

Members

NCM_NNLS_UMETHOD_NORMAL

Solve using normal equations and Cholesky decomposition

 

NCM_NNLS_UMETHOD_NORMAL_LU

Solve using normal equations and LU decomposition

 

NCM_NNLS_UMETHOD_QR

Solve using QR decomposition

 

NCM_NNLS_UMETHOD_DGELSD

Solve using QR decomposition (Lapack's dgelsd)

 

NCM_NNLS_UMETHOD_GSL

Solve using GSL's gsl_multifit_linear

 

NcmNNLS

typedef struct _NcmNNLS NcmNNLS;

Property Details

The “ncols” property

  “ncols”                    guint

Number of cols.

Owner: NcmNNLS

Flags: Read / Write / Construct Only

Allowed values: >= 1

Default value: 1


The “nrows” property

  “nrows”                    guint

Number of rows.

Owner: NcmNNLS

Flags: Read / Write / Construct Only

Allowed values: >= 1

Default value: 1


The “reltol” property

  “reltol”                   double

Relative tolerance.

Owner: NcmNNLS

Flags: Read / Write / Construct

Allowed values: [G_MINDOUBLE,0.1]

Default value: 2.22045e-16


The “umethod” property

  “umethod”                  NcmNNLSUMethod

Unconstrained method.

Owner: NcmNNLS

Flags: Read / Write / Construct

Default value: NCM_NNLS_UMETHOD_NORMAL