Linear solver in C++

Asked by Phil Marinier on 2010-05-28

I am attempting to use LinearSolver in C++ to do some simple matrix and vector solution, as in 3x3 or 5x5. I get this error.

[0]PETSC ERROR: --------------------- Error Message ------------------------------------
[0]PETSC ERROR: Object is in wrong state!
[0]PETSC ERROR: Not for unassembled matrix!
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Petsc Release Version 3.0.0, Patch 8, Fri Aug 21 14:02:12 CDT 2009

This seems to suggest that I cannot use LinearSolver in the code directly, that instead I have to make a form file. Is this in fact the case? and if so how would I make that form file?

This is my code in case it is relevant. I am doing some simple interpolation to calculate a step size for an iterative non-linear problem. This is the test program I wrote.

int main()
{
 //declare and initialise
 int degree = 4;
 double *x = new double[degree + 1];
 double *y = new double[degree + 1];

 std::vector<dolfin::uint> columns;
 std::vector<double> values;

 dolfin::LinearSolver linSolver;
 dolfin::Matrix X(degree + 1, degree + 1);
 dolfin::Vector Y(degree + 1);
 dolfin::Vector A(degree + 1);

 //generate test points
 x[0] = 0;
 y[0] = 0;

 x[1] = 1;
 y[1] = -2;

 x[2] = 2;
 y[2] = -1;

 x[3] = 3;
 y[3] = -2;

 x[4] = 4;
 y[4] = 0;

 //set vector
 Y.set_local(y);

 //set matrix
 for (int j = 0; j <= degree; j ++)
  columns.push_back(dolfin::uint(j));

 for (dolfin::uint row = 0; row < dolfin::uint(degree); row++){

  for (int j = 0; j <= degree; j++)
   values.push_back( power(x[row], (degree - j) ) );

  X.setrow(row, columns, values);

  values.clear();
 }

 //solve [X]{A} = {Y} where A is the coefficients of the polynomial
 linSolver.solve(X, A, Y);

 //output
 for (int j = 0; j <= degree; j++)
  std::cout << A[j] << std::endl;

 //clean up
 delete [] x;
 delete [] y;
}

Any help is appreciated. Thank you.

Question information

Language:
English Edit question
Status:
Solved
For:
DOLFIN Edit question
Assignee:
No assignee Edit question
Solved by:
Garth Wells
Solved:
2010-05-28
Last query:
2010-05-28
Last reply:
2010-05-28
Best Garth Wells (garth-wells) said : #1

On May 28 2010, Phil Marinier wrote:

>New question #112556 on DOLFIN:
>https://answers.launchpad.net/dolfin/+question/112556
>
> I am attempting to use LinearSolver in C++ to do some simple matrix and
> vector solution, as in 3x3 or 5x5. I get this error.
>
>
> [0]PETSC ERROR: --------------------- Error Message
> ------------------------------------ [0]PETSC ERROR: Object is in wrong
> state! [0]PETSC ERROR: Not for unassembled matrix! [0]PETSC ERROR:
> ------------------------------------------------------------------------
> [0]PETSC ERROR: Petsc Release Version 3.0.0, Patch 8, Fri Aug 21 14:02:12
> CDT 2009
>
>
> This seems to suggest that I cannot use LinearSolver in the code
> directly, that instead I have to make a form file. Is this in fact the
> case? and if so how would I make that form file?
>

You need to call A.apply(), X.apply(), etc before using the matrix and
vectors.

Manipulating small matrices using sparse linear algebra objects will be
dead slow. I would recommend using Armadillo (take a look in
dolfin/adaptivity for an example).

Garth

>
> This is my code in case it is relevant. I am doing some simple
> interpolation to calculate a step size for an iterative non-linear
> problem. This is the test program I wrote.
>
>
>int main()
>{
> //declare and initialise
> int degree = 4;
> double *x = new double[degree + 1];
> double *y = new double[degree + 1];
>
> std::vector<dolfin::uint> columns;
> std::vector<double> values;
>
> dolfin::LinearSolver linSolver;
> dolfin::Matrix X(degree + 1, degree + 1);
> dolfin::Vector Y(degree + 1);
> dolfin::Vector A(degree + 1);
>
> //generate test points
> x[0] = 0;
> y[0] = 0;
>
> x[1] = 1;
> y[1] = -2;
>
> x[2] = 2;
> y[2] = -1;
>
> x[3] = 3;
> y[3] = -2;
>
> x[4] = 4;
> y[4] = 0;
>
>
> //set vector
> Y.set_local(y);
>
> //set matrix
> for (int j = 0; j <= degree; j ++)
> columns.push_back(dolfin::uint(j));
>
> for (dolfin::uint row = 0; row < dolfin::uint(degree); row++){
>
> for (int j = 0; j <= degree; j++)
> values.push_back( power(x[row], (degree - j) ) );
>
> X.setrow(row, columns, values);
>
> values.clear();
> }
>
> //solve [X]{A} = {Y} where A is the coefficients of the polynomial
> linSolver.solve(X, A, Y);
>
> //output
> for (int j = 0; j <= degree; j++)
> std::cout << A[j] << std::endl;
>
> //clean up
> delete [] x;
> delete [] y;
>}
>
>
>Any help is appreciated. Thank you.
>
>You received this question notification because you are a member of
>DOLFIN Team, which is an answer contact for DOLFIN.
>
>_______________________________________________
>Mailing list: https://launchpad.net/~dolfin
>Post to : <email address hidden>
>Unsubscribe : https://launchpad.net/~dolfin
>More help : https://help.launchpad.net/ListHelp
>

Phil Marinier (lonewolf-13p) said : #2

Thanks Garth Wells, that solved my question.

Phil Marinier (lonewolf-13p) said : #3

Thanks. I will also check out Armadillo.

Nguyen Van Dang (dang-1032170) said : #4

Hello,
When running this program on Windows, I got the following error:
no matching function for call to 'dolfin:Vector::set_local(double*&)'
candidate is: virtual void dolfin::Vector::set_local(const dolfin::Array<double>&)
Can you help me to solve this problem?
Thanks in advance.
Nguyen Van Dang

Johan Hake (johan-hake) said : #5

On Tuesday May 31 2011 09:31:05 Nguyen Van Dang wrote:
> Question #112556 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/112556
>
> Nguyen Van Dang posted a new comment:
> Hello,
> When running this program on Windows, I got the following error:
> no matching function for call to 'dolfin:Vector::set_local(double*&)'
> candidate is: virtual void dolfin::Vector::set_local(const
> dolfin::Array<double>&) Can you help me to solve this problem?
> Thanks in advance.
> Nguyen Van Dang

You need to let x and y to be of type Array<double> instead of C-arrays.
Exchange these lines

 double *x = new double[degree + 1];
 double *y = new double[degree + 1];

with

 Array<double> x(degree + 1);
 Array<double> y(degree + 1);

and you can remove the cleaning of the arrays as well.

Johan

Nguyen Van Dang (dang-1032170) said : #6

I tried a simple program as following:
#include <iostream>
#include <dolfin.h>
int main()
{
 //declare and initialise
 int degree = 4;
 Array<double> x(degree + 1);
 Array<double> y(degree + 1);
 dolfin::Vector Y(degree + 1);
 //generate test points
 x[0] = 0;
 y[0] = 0;
 x[1] = 1;
 y[1] = -2;
 x[2] = 2;
 y[2] = -1;
 x[3] = 3;
 y[3] = -2;
 x[4] = 4;
 y[4] = 0;
 //set vector
 Y.set_local(y);
}
I got the error: 'Array' was not declared in this scope. Can you help me
again?
Thanks a lot.
Nguyen Van Dang

On Tue, May 31, 2011 at 6:55 PM, Johan Hake <
<email address hidden>> wrote:

> Question #112556 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/112556
>
> Johan Hake posted a new comment:
> On Tuesday May 31 2011 09:31:05 Nguyen Van Dang wrote:
> > Question #112556 on DOLFIN changed:
> > https://answers.launchpad.net/dolfin/+question/112556
> >
> > Nguyen Van Dang posted a new comment:
> > Hello,
> > When running this program on Windows, I got the following error:
> > no matching function for call to 'dolfin:Vector::set_local(double*&)'
> > candidate is: virtual void dolfin::Vector::set_local(const
> > dolfin::Array<double>&) Can you help me to solve this problem?
> > Thanks in advance.
> > Nguyen Van Dang
>
> You need to let x and y to be of type Array<double> instead of C-arrays.
> Exchange these lines
>
> double *x = new double[degree + 1];
> double *y = new double[degree + 1];
>
> with
>
> Array<double> x(degree + 1);
> Array<double> y(degree + 1);
>
> and you can remove the cleaning of the arrays as well.
>
> Johan
>
> --
> You received this question notification because you are a direct
> subscriber of the question.
>

--
Nguyen Van Dang

Home Tel: +84 7 92 21 36 28
Another email: <email address hidden>

Johan Hake (johan-hake) said : #7

Array is included in the dolfin namespace. Exchange Array with dolfin::Array,
or remove all dolfin:: from your code by including

  using namespace dolfin;

Have a look in any of the included C++ demos and you get a hint.

Johan

On Tuesday May 31 2011 10:15:57 Nguyen Van Dang wrote:
> Question #112556 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/112556
>
> Nguyen Van Dang posted a new comment:
> I tried a simple program as following:
> #include <iostream>
> #include <dolfin.h>
> int main()
> {
> //declare and initialise
> int degree = 4;
> Array<double> x(degree + 1);
> Array<double> y(degree + 1);
> dolfin::Vector Y(degree + 1);
> //generate test points
> x[0] = 0;
> y[0] = 0;
> x[1] = 1;
> y[1] = -2;
> x[2] = 2;
> y[2] = -1;
> x[3] = 3;
> y[3] = -2;
> x[4] = 4;
> y[4] = 0;
> //set vector
> Y.set_local(y);
> }
> I got the error: 'Array' was not declared in this scope. Can you help me
> again?
> Thanks a lot.
> Nguyen Van Dang
>
> On Tue, May 31, 2011 at 6:55 PM, Johan Hake <
>
> <email address hidden>> wrote:
> > Question #112556 on DOLFIN changed:
> > https://answers.launchpad.net/dolfin/+question/112556
> >
> > Johan Hake posted a new comment:
> >
> > On Tuesday May 31 2011 09:31:05 Nguyen Van Dang wrote:
> > > Question #112556 on DOLFIN changed:
> > > https://answers.launchpad.net/dolfin/+question/112556
> > >
> > > Nguyen Van Dang posted a new comment:
> > > Hello,
> > > When running this program on Windows, I got the following error:
> > > no matching function for call to 'dolfin:Vector::set_local(double*&)'
> > > candidate is: virtual void dolfin::Vector::set_local(const
> > > dolfin::Array<double>&) Can you help me to solve this problem?
> > > Thanks in advance.
> > > Nguyen Van Dang
> >
> > You need to let x and y to be of type Array<double> instead of C-arrays.
> > Exchange these lines
> >
> > double *x = new double[degree + 1];
> > double *y = new double[degree + 1];
> >
> > with
> >
> > Array<double> x(degree + 1);
> > Array<double> y(degree + 1);
> >
> > and you can remove the cleaning of the arrays as well.
> >
> > Johan
> >
> > --
> > You received this question notification because you are a direct
> > subscriber of the question.

Nguyen Van Dang (dang-1032170) said : #8

Thanks Mr. Johan Hake. That solved my problem.

On Tue, May 31, 2011 at 7:31 PM, Johan Hake <
<email address hidden>> wrote:

> Question #112556 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/112556
>
> Johan Hake posted a new comment:
> Array is included in the dolfin namespace. Exchange Array with
> dolfin::Array,
> or remove all dolfin:: from your code by including
>
> using namespace dolfin;
>
> Have a look in any of the included C++ demos and you get a hint.
>
> Johan
>
> On Tuesday May 31 2011 10:15:57 Nguyen Van Dang wrote:
> > Question #112556 on DOLFIN changed:
> > https://answers.launchpad.net/dolfin/+question/112556
> >
> > Nguyen Van Dang posted a new comment:
> > I tried a simple program as following:
> > #include <iostream>
> > #include <dolfin.h>
> > int main()
> > {
> > //declare and initialise
> > int degree = 4;
> > Array<double> x(degree + 1);
> > Array<double> y(degree + 1);
> > dolfin::Vector Y(degree + 1);
> > //generate test points
> > x[0] = 0;
> > y[0] = 0;
> > x[1] = 1;
> > y[1] = -2;
> > x[2] = 2;
> > y[2] = -1;
> > x[3] = 3;
> > y[3] = -2;
> > x[4] = 4;
> > y[4] = 0;
> > //set vector
> > Y.set_local(y);
> > }
> > I got the error: 'Array' was not declared in this scope. Can you help me
> > again?
> > Thanks a lot.
> > Nguyen Van Dang
> >
> > On Tue, May 31, 2011 at 6:55 PM, Johan Hake <
> >
> > <email address hidden>> wrote:
> > > Question #112556 on DOLFIN changed:
> > > https://answers.launchpad.net/dolfin/+question/112556
> > >
> > > Johan Hake posted a new comment:
> > >
> > > On Tuesday May 31 2011 09:31:05 Nguyen Van Dang wrote:
> > > > Question #112556 on DOLFIN changed:
> > > > https://answers.launchpad.net/dolfin/+question/112556
> > > >
> > > > Nguyen Van Dang posted a new comment:
> > > > Hello,
> > > > When running this program on Windows, I got the following error:
> > > > no matching function for call to 'dolfin:Vector::set_local(double*&)'
> > > > candidate is: virtual void dolfin::Vector::set_local(const
> > > > dolfin::Array<double>&) Can you help me to solve this problem?
> > > > Thanks in advance.
> > > > Nguyen Van Dang
> > >
> > > You need to let x and y to be of type Array<double> instead of
> C-arrays.
> > > Exchange these lines
> > >
> > > double *x = new double[degree + 1];
> > > double *y = new double[degree + 1];
> > >
> > > with
> > >
> > > Array<double> x(degree + 1);
> > > Array<double> y(degree + 1);
> > >
> > > and you can remove the cleaning of the arrays as well.
> > >
> > > Johan
> > >
> > > --
> > > You received this question notification because you are a direct
> > > subscriber of the question.
>
> --
> You received this question notification because you are a direct
> subscriber of the question.
>

--
Nguyen Van Dang

Home Tel: +84 7 92 21 36 28
Another email: <email address hidden>