I have looked into the problem ore and determined this is more likely a problem with how "ydot" is being handled. I checked my assignment to ydot within the RHS function and the contents are correct. However, after adding a print statement in the CVODE code after the call to
retval = f(tn, zn[0], zn[1], user_data);
The values of "ydot" or "zn[1]" are much different than those that should to be output (based off what I am setting in the RHS function).
Based on what I've read from the documentation it seems that my definition of ydot in the RHS function should be adequate but perhaps I am missing something.
Thanks,
Brian
bowens1 wrote:
Hello,
I am currently implementing CVODE into an in-house FEA code written in C++
that makes thorough use of object oriented programming. Within a particular
method of a class I am trying to implement CVODE. I have made my RHS
function a static function outside the class as I have seen done in the
examples. I am having some problem in that at the second call to the
RHSfunction the values passed in under "N_vector y" seem clearly incorrect.
For example my IVP for the problem is [0.79 0 0] and the second call to RHS
function has the y vector as [-5.0784477968231915e+057
-5.0784477968231915e+057 -5.0784477968231915e+057]. It seems as though the y
vector is being corrupted some where in the main CVODE function.
To test out sundials I wrote a smaller matlab code that is not nearly as
robust as the C++ implementation. The matlab implementation worked with
great success. I am using the same settings in both C++ and matlab
implementations. Also the two implementations are giving identical output
during the first call to the RHS function, the C++ code gives erroneous
results after this.
If needed here is my RHS function:
DiffusionModel is a class and user_data points to an object of this class
that I later assign to dmod. The method "solve_for_qdot()" solves for qdot
and assigns it to some data object "udot" local to the DiffusionModel class.
//==============================================================
static int RHSfunction(realtype t, N_Vector y, N_Vector ydot, void
*user_data)
{
DiffusionModel* dmod;
dmod =(DiffusionModel*)user_data;
double *qdot;
double *q;
q = NV_DATA_S(y);
dmod->solve_for_qdot(t,q);
qdot = dmod->udot;
ydot = N_VMake_Serial(dmod->totalNumDof,qdot);
int flag=0;
return flag;
}
//==============================================================
I greatly appreciate any help in this matter, and I'll be glad to provide
anymore details as needed.
Thanks,
Brian
