Author |
Topic |
|
ggarrett
Canada
3 Posts |
Posted - 03/17/2015 : 11:28:05 AM
|
Origin Ver. and Service Release (Select Help-->About Origin): OriginPro 9.0.0 (64-bit) SR2 b87 Operating System: Windows 7 64 bit
I'm trying to modify the ODE fitting tutorial found here (http://www.originlab.com/doc/Tutorials/Fitting-Ordinary-Differential-Equation) to something a bit more useful. I went through the tutorial and got the original programming working fine.
The ODE I'm trying to fit is below. I'm trying to fit concentration vs time data. D, CatT, and y0 are fixed parameters. I'm trying to fit k1 and k2.
dy(t)/dt = (k1/D)*(D*y[0]+y[0]^2)/(1.0+(k1+k2)/(k2*D))*CatT;
I've attached the code below. It compiles fine, but doesn't evaluate. A previous version I had of it did evaluate successfully but I get a similar result. When I run the code I get no fit to the data, y remains fixed at y0. If anyone can give me a hand I would be grateful.
#pragma numlittype(push, TRUE) #pragma numlittype(push, TRUE) #pragma numlittype(push, TRUE) #pragma numlittype(push, TRUE) #pragma numlittype(push, TRUE) #pragma numlittype(push, TRUE) #pragma numlittype(push, TRUE) #pragma warning(error : 15618) #include <origin.h> #include <oc_Nag8.h> #include <ONLSF.H>
struct user // parameter in the ODE { double k1; double k2; double D; double CatT; }; //Define the differentiate equation: y'=a*y static void NAG_CALL f(Integer neq, double t, double y[], double yp[], Nag_User *comm) { neq; //Number of ordinary differential equations t; //Independent variable y; //Dependent variables yp; //First derivative struct user *sp = (struct user *)(comm->p); double k1; double k2; double D; double CatT; k1 = sp->k1; k2 = sp->k2; D = sp->D; CatT = sp->CatT; yp[0] = (k1/D)*(D*y[0]+y[0]^2)/(1.0+(k1+k2)/(k2*D))*CatT; } //Use Runge-Kutta ODE23 to solve ODE static bool nag_ode_fit( const double k1, const double k2, const double D, const double CatT, const double y0, const double tstart, const double tend, const int nout, vector &vP ) { //nout: Number of points to output if( nout < 2 ) return false; vP.SetSize( nout ); vP[0] = y0; int neq = 1; //Number of ordinary differential equations Nag_RK_method method; double hstart, tgot, tinc; double tol, twant; int i, j; vector thres(neq), ygot(neq), ymax(neq), ypgot(neq), ystart(neq); Nag_ErrorAssess errass; Nag_ODE_RK opt; Nag_User comm; struct user s; s.k1 = k1; s.k2 = k2; s.D = D; s.CatT = CatT; comm.p = (Pointer)&s; ystart[0] = y0; for (i=0; i<neq; i++) thres[i] = 1.0e-8; errass = Nag_ErrorAssess_off; hstart = 0; method = Nag_RK_2_3; tinc = (tend-tstart)/(nout-1); tol = 1.0e-3; NagError nagErr1; //Setup ODE d02pvc(neq, tstart, ystart, tend, tol, thres, method, Nag_RK_range, errass, hstart, &opt, &nagErr1); if( nagErr1.code != NE_NOERROR ) return false; for (j=1; j<nout; j++) { twant = tstart + j*tinc; NagError nagErr2; //Solve ODE d02pcc(neq, f, twant, &tgot, ygot, ypgot, ymax, &opt, &comm, &nagErr2); if( nagErr2.code != NE_NOERROR ) return false; vP[j] = ygot[0]; } //Free functions for Runge–Kutta suite d02ppc(&opt); return true; }
//---------------------------------------------------------- // void _nlsfBlackmondFit( // Fit Parameter(s): double k1, double k2, double D, double CatT, double y0, // Independent Variable(s): double t, // Dependent Variable(s): double& y) { // Beginning of editable part NLFitContext *pCtxt = Project.GetNLFitContext(); if ( pCtxt ) { static vector vX, vY; static int nSize; BOOL bIsNewParamValues = pCtxt->IsNewParamValues(); // If parameters were updated, we will recalculate the ODE result. if ( bIsNewParamValues ) { //Initial and final values of the independent variable double tstart = 0, tend = 69.3, tinc; int nout = 300; //Number of points tinc = (tend-tstart)/(nout-1); vX.Data( tstart, tend, tinc ); nSize = vX.GetSize(); if( !nag_ode_fit( k1, k2, D, CatT, y0, tstart, tend, nout, vY) ) return; } //Interpolate y from fitting data's x on the ODE result. ocmath_interpolate( &t, &y, 1, vX, vY, nSize ); } // End of editable part } |
|
ggarrett
Canada
3 Posts |
Posted - 03/17/2015 : 12:33:00 PM
|
I lied, code works. Put silly restrictions onto my parameters |
|
|
ggarrett
Canada
3 Posts |
Posted - 03/17/2015 : 12:36:07 PM
|
Except I still don't get a result when I hit evaluate. S[its out y = -- |
|
|
cc261
22 Posts |
Posted - 03/20/2015 : 05:34:02 AM
|
for ODE fitting problem, try a software named 1stOpt or Auto2Fit, 100-times easy than Origin |
|
|
|
Topic |
|
|
|