Upgraded from 1.2 to 1.4, now I get a PETSc error on 'VecCopy'

Question

I updated from 1.2 to 1.4 and now on my call to

*x.get()->vector() = *other.x.get()->vector()->copy();

I get an error from PETSC:

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
***
***     fenics@fenicsproject.org
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error:   Unable to successfully call PETSc function 'VecCopy'.
*** Reason:  PETSc error code is: 73.
*** Where:   This error was encountered inside ~/Work/FEniCS/src/dolfin-1.4.0/dolfin/la/PETScVector.cpp.
*** Process: unknown
*** 
*** DOLFIN version: 1.4.0
*** Git changeset:  
*** -------------------------------------------------------------------------

The PETSc error

[0]PETSC ERROR: --------------------- Error Message ------------------------------------
[0]PETSC ERROR: Object is in wrong state!
[0]PETSC ERROR: Not for unassembled vector!
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Petsc Release Version 3.4.2, Jul, 02, 2013 
[0]PETSC ERROR: See docs/changes/index.html for recent updates.
[0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
[0]PETSC ERROR: See docs/index.html for manual pages.
[0]PETSC ERROR: ------------------------------------------------------------------------

Do I need to do something differently with 1.4? Perhaps an apply call?

---

This seems to work, though rather syntactically heavy:

std::vector<double> temp_val;
other.x->vector()->get_local(temp_val);
x->vector()->set_local(temp_val);
std::dynamic_pointer_cast<PETScVector>(x->vector())->apply("");

It would also be nice to not need the working vector temp_val to perform the copy

Answer 1

The solution I found was to do the following

std::vector<double> temp_val;
std::dynamic_pointer_cast<PETScVector>(other.x->vector())->apply("");
other.x->vector()->get_local(temp_val);
x->vector()->set_local(temp_val);
std::dynamic_pointer_cast<PETScVector>(x->vector())->apply("");

---

Note the dynamic_pointer_casts are not needed, vector() itself has apply("").

Comments

Hi!

First of all thank you for the solution! It works!

Could you explain what does a line of code below mean?

std::dynamic_pointer_cast<PETScVector>(x->vector())->apply("");

Why we have to put the line of code before

other.x->vector()->get_local(temp_val);

and after

x->vector()->set_local(temp_val);

Thanks in advance!

---

Hey Maks,

Because the type of vector being used is not known (it's type is just 'some' generic vector) we use the fact that we know it happens to be a PETScVector to invoke a PETScVector method "apply". If it wasn't a PETScVector this would fail badly at runtime.

The apply forces the communication to occur so that the local in memory values are 'up to date'. With updated local values in memory we do the copy locally, and then we again call the apply to communicate the changes.

---

Hi Charles!

Thanks for the answer!

I am partially get it.

Could we call the apply("") without cast to a PETScVector as a virtual function?

I understand that we have to call apply("") after set_local(temp_val) to finalize assembly of a tensor.

Why we should call it before get_local()?
We doesn't change a tensor by this call.

---

It's been a few releases since I dove into this, so it might be possible to just call apply without the cast (I think it wasnt available in general at the time).

If you don't call it before get_local(), you may not have the data from the remote processes. If you are running without MPI you should be able to skip it.

---

Yes, It works without cast to the PETScVector. I have checked it.

Thanks for notice about MPI!

In case of MPI it gets all data on a local process, so if I will have 2 processes then I have to call get_local() and set_local() only on main process or on each one and every process get it's peace of data?

---

With MPI you each process should be doing the same thing. So the get_local and set_local would need to be done on each process. It will handle the communication for you, thankfully.

---

Thank you very much!