00001 #include<iostream>
00002 #include<fstream>
00003 #include<math.h>
00004 #include<string.h>
00005 #include<stdio.h>
00006
00007 #include "../def.h"
00008 #include "../model.h"
00009 #include "../nr.h"
00010
00011 using namespace std;
00012
00013
00014
00015
00016
00017 void tabulateValues (Glob *globs,GlobExp *ex,double t0, double t1, double *t,
00018 long n, double *state, double **y, double ***dmds,
00019 double ***dmdp, double **dyds,double **dydp);
00020
00021 void simInit(GlobExp *ex,Glob *globs)
00022 {
00023 long i,j;
00024 long nvar=ex->nvar;
00025 long nPoints=ex->nPoints;
00026 long idum=-42;
00027 double *t=dvector(1,2);
00028 double **y=dmatrix(1,2,1,nvar);
00029
00030 if (ex->y0)
00031 {
00032 for (i=1; i<=nvar; ++i)
00033 ex->yTry[1][i]=ex->y0[i];
00034 }
00035 else
00036 {
00037 for (i=1; i<=nvar; ++i)
00038 {
00039 ex->yTry[1][i]=DefYValues[i-1];
00040 }
00041 }
00042
00043
00044
00045 for(i=2;i<=nPoints;i++)
00046 {
00047 for(j=1;j<=nvar;j++)
00048 ex->yTry[i][j]=ex->yTry[i-1][j];
00049 t[1]=ex->mesh[i-1];
00050 t[2]=ex->mesh[i];
00051 tabulateValues(globs,ex,t[1],t[2],t,2,ex->yTry[i],y,NULL,NULL,NULL,NULL);
00052 for(j=1;j<=nvar;j++)
00053 ex->yTry[i][j]+=globs->pert*ex->yTry[i][j]*gasdev(&idum);
00054 }
00055
00056
00057
00058 free_dvector(t,1,2);
00059 free_dmatrix(y,1,2,1,nvar);
00060 }