TPCCLIB
Loading...
Searching...
No Matches
simulate.c
Go to the documentation of this file.
1
7/*****************************************************************************/
8#include "libtpcmodel.h"
9/*****************************************************************************/
10
11/*****************************************************************************/
29 double *t,
31 double *ca,
33 int nr,
35 double k1,
37 double k2,
39 double k3,
41 double k4,
43 double k5,
45 double k6,
47 double *ct,
49 double *cta,
51 double *ctb,
53 double *ctc
54) {
55 double b, c, d, w, z, dt2;
56 double cai, ca_last, t_last;
57 double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
58 double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
59
60
61 /* Check for data */
62 if(nr<2) return 1;
63 if(ct==NULL) return 2;
64
65 /* Check actual parameter number */
66 if(k1<0.0) return 3;
67 if(k3<=0.0) {k3=0.0; if(k2<=0.0) {k2=0.0;}}
68 else if(k5<=0.0) {k5=0.0; if(k4<=0.0) {k4=0.0;}}
69 else if(k6<=0.0) {k6=0.0;}
70
71 /* Calculate curves */
72 t_last=0.0; if(t[0]<t_last) t_last=t[0];
73 cai=ca_last=0.0;
74 ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=0.0;
75 ct1=ct2=ct3=ct1i=ct2i=ct3i=0.0;
76 for(int i=0; i<nr; i++) {
77 /* delta time / 2 */
78 dt2=0.5*(t[i]-t_last);
79 /* calculate values */
80 if(dt2<0.0) {
81 return 5;
82 } else if(dt2>0.0) {
83 /* arterial integral */
84 cai+=(ca[i]+ca_last)*dt2;
85 /* partial results */
86 b=ct1i_last+dt2*ct1_last;
87 c=ct2i_last+dt2*ct2_last;
88 d=ct3i_last+dt2*ct3_last;
89 w=k4 + k5 - (k5*k6*dt2)/(1.0+k6*dt2);
90 z=1.0+w*dt2;
91 /* 1st tissue compartment and its integral */
92 ct1 = (
93 + k1*z*cai + (k3*k4*dt2 - (k2+k3)*z)*b
94 + k4*c + k4*k6*dt2*d/(1.0+k6*dt2)
95 ) / ( z*(1.0 + dt2*(k2+k3)) - k3*k4*dt2*dt2 );
96 ct1i = ct1i_last + dt2*(ct1_last+ct1);
97 /* 2nd tissue compartment and its integral */
98 ct2 = (k3*ct1i - w*c + k6*d/(1.0+k6*dt2)) / z;
99 ct2i = ct2i_last + dt2*(ct2_last+ct2);
100 /* 3rd tissue compartment and its integral */
101 ct3 = (k5*ct2i - k6*d) / (1.0 + k6*dt2);
102 ct3i = ct3i_last + dt2*(ct3_last+ct3);
103 }
104 /* copy values to argument arrays; set very small values to zero */
105 ct[i]=ct1+ct2+ct3; if(fabs(ct[i])<1.0e-12) ct[i]=0.0;
106 if(cta!=NULL) {cta[i]=ct1; if(fabs(cta[i])<1.0e-12) cta[i]=0.0;}
107 if(ctb!=NULL) {ctb[i]=ct2; if(fabs(ctb[i])<1.0e-12) ctb[i]=0.0;}
108 if(ctc!=NULL) {ctc[i]=ct3; if(fabs(ctc[i])<1.0e-12) ctc[i]=0.0;}
109 /* prepare to the next loop */
110 t_last=t[i]; ca_last=ca[i];
111 ct1_last=ct1; ct1i_last=ct1i;
112 ct2_last=ct2; ct2i_last=ct2i;
113 ct3_last=ct3; ct3i_last=ct3i;
114 }
115
116 return 0;
117}
118/*****************************************************************************/
119
120/*****************************************************************************/
138 double *t,
140 double *ca,
142 int nr,
144 double k1,
146 double k2,
148 double k3,
150 double k4,
152 double k5,
154 double k6,
156 double *ct,
158 double *cta,
160 double *ctb,
162 double *ctc
163) {
164 int i;
165 double dt2, r, s, u, v, w;
166 double cai, ca_last, t_last;
167 double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
168 double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
169
170
171 /* Check for data */
172 if(nr<2) return 1;
173 if(ct==NULL) return 2;
174
175 /* Check parameters */
176 if(k1<0.0) return 3;
177
178 /* Calculate curves */
179 t_last=0.0; if(t[0]<t_last) t_last=t[0];
180 cai=ca_last=0.0;
181 ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=0.0;
182 ct1=ct2=ct3=ct1i=ct2i=ct3i=0.0;
183 for(i=0; i<nr; i++) {
184 /* delta time / 2 */
185 dt2=0.5*(t[i]-t_last);
186 /* calculate values */
187 if(dt2<0.0) {
188 return 5;
189 } else if(dt2>0.0) {
190 /* arterial integral */
191 cai+=(ca[i]+ca_last)*dt2;
192 /* Calculate partial results */
193 r=1.0+k4*dt2;
194 s=1.0+k6*dt2;
195 u=ct1i_last+dt2*ct1_last;
196 v=ct2i_last+dt2*ct2_last;
197 w=ct3i_last+dt2*ct3_last;
198 /* 1st tissue compartment and its integral */
199 ct1 = ( k1*cai - (k2 + (k3/r) + (k5/s))*u + (k4/r)*v + (k6/s)*w )
200 / ( 1.0 + dt2*(k2 + (k3/r) + (k5/s)) );
201 ct1i = ct1i_last + dt2*(ct1_last+ct1);
202 /* 2nd tissue compartment and its integral */
203 ct2 = (k3*ct1i - k4*v) / r;
204 ct2i = ct2i_last + dt2*(ct2_last+ct2);
205 /* 3rd tissue compartment and its integral */
206 ct3 = (k5*ct1i - k6*w) / s;
207 ct3i = ct3i_last + dt2*(ct3_last+ct3);
208 }
209 /* copy values to argument arrays; set very small values to zero */
210 ct[i]=ct1+ct2+ct3; if(fabs(ct[i])<1.0e-12) ct[i]=0.0;
211 if(cta!=NULL) {cta[i]=ct1; if(fabs(cta[i])<1.0e-12) cta[i]=0.0;}
212 if(ctb!=NULL) {ctb[i]=ct2; if(fabs(ctb[i])<1.0e-12) ctb[i]=0.0;}
213 if(ctc!=NULL) {ctc[i]=ct3; if(fabs(ctc[i])<1.0e-12) ctc[i]=0.0;}
214 /* prepare to the next loop */
215 t_last=t[i]; ca_last=ca[i];
216 ct1_last=ct1; ct1i_last=ct1i;
217 ct2_last=ct2; ct2i_last=ct2i;
218 ct3_last=ct3; ct3i_last=ct3i;
219 }
220
221 return 0;
222}
223/*****************************************************************************/
224
225/*****************************************************************************/
245 double *t,
247 double *ca,
249 double *cb,
251 int nr,
253 double k1,
255 double k2,
257 double k3,
259 double k4,
261 double k5,
263 double k6,
265 double f,
267 double vb,
269 double fa,
271 double *cpet,
273 double *cta,
275 double *ctb,
277 double *ctc,
279 double *ctab,
281 double *ctvb
282) {
283 int i;
284 double b, c, d, w, z, dt2, va, vv;
285 double cai, ca_last, t_last, dct, cvb;
286 double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
287 double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
288
289
290 /* Check for data */
291 if(nr<2) return 1;
292 if(cpet==NULL) return 2;
293
294 /* Check parameters */
295 if(k1<0.0) return 3;
296 if(vb<0.0 || vb>=1.0) return 4;
297 if(fa<=0.0 || fa>1.0) return 5;
298 va=fa*vb; vv=(1.0-fa)*vb;
299
300 /* Calculate curves */
301 t_last=0.0; if(t[0]<t_last) t_last=t[0];
302 cai=ca_last=0.0;
303 ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=0.0;
304 ct1=ct2=ct3=ct1i=ct2i=ct3i=0.0;
305 for(i=0; i<nr; i++) {
306 /* delta time / 2 */
307 dt2=0.5*(t[i]-t_last);
308 /* calculate values */
309 if(dt2<0.0) {
310 return 5;
311 } else if(dt2>0.0) {
312 /* arterial integral */
313 cai+=(ca[i]+ca_last)*dt2;
314 /* partial results */
315 b=ct1i_last+dt2*ct1_last;
316 c=ct2i_last+dt2*ct2_last;
317 d=ct3i_last+dt2*ct3_last;
318 w=k4 + k5 - (k5*k6*dt2)/(1.0+k6*dt2);
319 z=1.0+w*dt2;
320 /* 1st tissue compartment and its integral */
321 ct1 = (
322 + k1*z*cai + (k3*k4*dt2 - (k2+k3)*z)*b
323 + k4*c + k4*k6*dt2*d/(1.0+k6*dt2)
324 ) / ( z*(1.0 + dt2*(k2+k3)) - k3*k4*dt2*dt2 );
325 ct1i = ct1i_last + dt2*(ct1_last+ct1);
326 /* 2nd tissue compartment and its integral */
327 ct2 = (k3*ct1i - w*c + k6*d/(1.0+k6*dt2)) / z;
328 ct2i = ct2i_last + dt2*(ct2_last+ct2);
329 /* 3rd tissue compartment and its integral */
330 ct3 = (k5*ct2i - k6*d) / (1.0 + k6*dt2);
331 ct3i = ct3i_last + dt2*(ct3_last+ct3);
332 }
333 /* Venous curve */
334 if(f>0.) {dct = k1*ca[i] - k2*ct1; cvb = cb[i] - dct/f;} else cvb=cb[i];
335 /* copy values to argument arrays; set very small values to zero */
336 cpet[i]= va*cb[i] + vv*cvb + (1.0-vb)*(ct1+ct2+ct3);
337 if(fabs(cpet[i])<1.0e-12) cpet[i]=0.0;
338 if(cta!=NULL) {cta[i]=(1.0-vb)*ct1; if(fabs(cta[i])<1.0e-12) cta[i]=0.0;}
339 if(ctb!=NULL) {ctb[i]=(1.0-vb)*ct2; if(fabs(ctb[i])<1.0e-12) ctb[i]=0.0;}
340 if(ctc!=NULL) {ctc[i]=(1.0-vb)*ct3; if(fabs(ctc[i])<1.0e-12) ctc[i]=0.0;}
341 if(ctab!=NULL) {ctab[i]=va*cb[i]; if(fabs(ctab[i])<1.0e-12) ctab[i]=0.0;}
342 if(ctvb!=NULL) {ctvb[i]=vv*cvb; if(fabs(ctvb[i])<1.0e-12) ctvb[i]=0.0;}
343 /* prepare to the next loop */
344 t_last=t[i]; ca_last=ca[i];
345 ct1_last=ct1; ct1i_last=ct1i;
346 ct2_last=ct2; ct2i_last=ct2i;
347 ct3_last=ct3; ct3i_last=ct3i;
348 }
349
350 return 0;
351}
352/*****************************************************************************/
353
354/*****************************************************************************/
375 double *t,
377 double *ca,
379 double *cb,
381 int nr,
383 double k1,
385 double k2,
387 double k3,
389 double k4,
391 double k5,
393 double k6,
395 double f,
397 double vb,
399 double fa,
401 double *cpet,
403 double *cta,
405 double *ctb,
407 double *ctc,
409 double *ctab,
411 double *ctvb
412) {
413 int i;
414 double dt2, r, s, u, v, w, va, vv;
415 double cai, ca_last, t_last, dct, cvb;
416 double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
417 double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
418
419
420 /* Check for data */
421 if(nr<2) return 1;
422 if(cpet==NULL) return 2;
423
424 /* Check parameters */
425 if(k1<0.0) return 3;
426 if(vb<0.0 || vb>=1.0) return 4;
427 if(fa<=0.0 || fa>1.0) return 5;
428 va=fa*vb; vv=(1.0-fa)*vb;
429
430 /* Calculate curves */
431 t_last=0.0; if(t[0]<t_last) t_last=t[0];
432 cai=ca_last=0.0;
433 ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=0.0;
434 ct1=ct2=ct3=ct1i=ct2i=ct3i=0.0;
435 for(i=0; i<nr; i++) {
436 /* delta time / 2 */
437 dt2=0.5*(t[i]-t_last);
438 /* calculate values */
439 if(dt2<0.0) {
440 return 5;
441 } else if(dt2>0.0) {
442 /* arterial integral */
443 cai+=(ca[i]+ca_last)*dt2;
444 /* Calculate partial results */
445 r=1.0+k4*dt2;
446 s=1.0+k6*dt2;
447 u=ct1i_last+dt2*ct1_last;
448 v=ct2i_last+dt2*ct2_last;
449 w=ct3i_last+dt2*ct3_last;
450 /* 1st tissue compartment and its integral */
451 ct1 = ( k1*cai - (k2 + (k3/r) + (k5/s))*u + (k4/r)*v + (k6/s)*w )
452 / ( 1.0 + dt2*(k2 + (k3/r) + (k5/s)) );
453 ct1i = ct1i_last + dt2*(ct1_last+ct1);
454 /* 2nd tissue compartment and its integral */
455 ct2 = (k3*ct1i - k4*v) / r;
456 ct2i = ct2i_last + dt2*(ct2_last+ct2);
457 /* 3rd tissue compartment and its integral */
458 ct3 = (k5*ct1i - k6*w) / s;
459 ct3i = ct3i_last + dt2*(ct3_last+ct3);
460 }
461 /* Venous curve */
462 if(f>0.) {dct = k1*ca[i] - k2*ct1; cvb = cb[i] - dct/f;} else cvb=cb[i];
463 /* copy values to argument arrays; set very small values to zero */
464 cpet[i]= va*cb[i] + vv*cvb + (1.0-vb)*(ct1+ct2+ct3);
465 if(fabs(cpet[i])<1.0e-12) cpet[i]=0.0;
466 if(cta!=NULL) {cta[i]=(1.0-vb)*ct1; if(fabs(cta[i])<1.0e-12) cta[i]=0.0;}
467 if(ctb!=NULL) {ctb[i]=(1.0-vb)*ct2; if(fabs(ctb[i])<1.0e-12) ctb[i]=0.0;}
468 if(ctc!=NULL) {ctc[i]=(1.0-vb)*ct3; if(fabs(ctc[i])<1.0e-12) ctc[i]=0.0;}
469 if(ctab!=NULL) {ctab[i]=va*cb[i]; if(fabs(ctab[i])<1.0e-12) ctab[i]=0.0;}
470 if(ctvb!=NULL) {ctvb[i]=vv*cvb; if(fabs(ctvb[i])<1.0e-12) ctvb[i]=0.0;}
471 /* prepare to the next loop */
472 t_last=t[i]; ca_last=ca[i];
473 ct1_last=ct1; ct1i_last=ct1i;
474 ct2_last=ct2; ct2i_last=ct2i;
475 ct3_last=ct3; ct3i_last=ct3i;
476 }
477
478 return 0;
479}
480/*****************************************************************************/
481
482/*****************************************************************************/
502 double *t,
504 double *ca,
506 int nr,
508 double k1,
510 double k2,
512 double k3,
514 double kLoss,
516 double *ct,
518 double *cta,
520 double *ctb
521) {
522 int i;
523 double b, c, dt2;
524 double cai, ca_last, t_last;
525 double ct1, ct1_last, ct2, ct2_last;
526 double ct1i, ct1i_last, ct2i, ct2i_last;
527
528
529 /* Check for data */
530 if(nr<2) return 1;
531 if(ct==NULL) return 2;
532
533 /* Check actual parameter number */
534 if(k1<0.0) return 3;
535
536 /* Calculate curves */
537 t_last=0.0; if(t[0]<t_last) t_last=t[0];
538 cai=ca_last=0.0;
539 ct1_last=ct2_last=ct1i_last=ct2i_last=ct1=ct2=ct1i=ct2i=0.0;
540 for(i=0; i<nr; i++) {
541 /* delta time / 2 */
542 dt2=0.5*(t[i]-t_last);
543 /* calculate values */
544 if(dt2<0.0) {
545 return 5;
546 } else if(dt2>0.0) {
547 /* arterial integral */
548 cai+=(ca[i]+ca_last)*dt2;
549 /* partial results */
550 b=ct1i_last+dt2*ct1_last;
551 c=ct2i_last+dt2*ct2_last;
552 /* 1st tissue compartment and its integral */
553 ct1 = (k1*cai - (k2+k3)*b) / (1.0 + (k2+k3)*dt2 );
554 ct1i = ct1i_last + dt2*(ct1_last+ct1);
555 /* 2nd tissue compartment and its integral */
556 ct2 = (k3*ct1i - kLoss*c) / (1.0 + kLoss*dt2);
557 ct2i = ct2i_last + dt2*(ct2_last+ct2);
558 }
559 /* copy values to argument arrays; set very small values to zero */
560 ct[i]=ct1+ct2; if(fabs(ct[i])<1.0e-12) ct[i]=0.0;
561 if(cta!=NULL) {cta[i]=ct1; if(fabs(cta[i])<1.0e-12) cta[i]=0.0;}
562 if(ctb!=NULL) {ctb[i]=ct2; if(fabs(ctb[i])<1.0e-12) ctb[i]=0.0;}
563 /* prepare to the next loop */
564 t_last=t[i]; ca_last=ca[i];
565 ct1_last=ct1; ct1i_last=ct1i;
566 ct2_last=ct2; ct2i_last=ct2i;
567 }
568
569 return 0;
570}
571/*****************************************************************************/
572
573/*****************************************************************************/
594 double *t,
596 double *ca,
598 double *cb,
600 int nr,
602 double k1,
604 double k2,
606 double k3,
608 double kL,
610 double f,
612 double vb,
614 double fa,
616 double *cpet,
618 double *cta,
620 double *ctb,
622 double *ctab,
624 double *ctvb
625) {
626 int i;
627 double dt2, b, c, va, vv;
628 double cai, ca_last, t_last, dct, cvb;
629 double ct1, ct1_last, ct2, ct2_last;
630 double ct1i, ct1i_last, ct2i, ct2i_last;
631
632
633 /* Check for data */
634 if(nr<2) return 1;
635 if(cpet==NULL) return 2;
636
637 /* Check parameters */
638 if(k1<0.0) return 3;
639 if(vb<0.0 || vb>=1.0) return 4;
640 if(fa<=0.0 || fa>1.0) return 5;
641 va=fa*vb; vv=(1.0-fa)*vb;
642
643 /* Calculate curves */
644 t_last=0.0; if(t[0]<t_last) t_last=t[0];
645 cai=ca_last=0.0;
646 ct1_last=ct2_last=ct1i_last=ct2i_last=ct1=ct2=ct1i=ct2i=0.0;
647 for(i=0; i<nr; i++) {
648 /* delta time / 2 */
649 dt2=0.5*(t[i]-t_last);
650 /* calculate values */
651 if(dt2<0.0) {
652 return 5;
653 } else if(dt2>0.0) {
654 /* arterial integral */
655 cai+=(ca[i]+ca_last)*dt2;
656 /* Calculate partial results */
657 b=ct1i_last+dt2*ct1_last;
658 c=ct2i_last+dt2*ct2_last;
659 /* 1st tissue compartment and its integral */
660 ct1 = (k1*cai - (k2+k3)*b) / (1.0 + (k2+k3)*dt2 );
661 ct1i = ct1i_last + dt2*(ct1_last+ct1);
662 /* 2nd tissue compartment and its integral */
663 ct2 = (k3*ct1i - kL*c) / (1.0 + kL*dt2);
664 ct2i = ct2i_last + dt2*(ct2_last+ct2);
665 }
666 /* Venous curve */
667 if(f>0.) {dct = k1*ca[i] - k2*ct1 - kL*ct2; cvb = cb[i] - dct/f;}
668 else cvb=cb[i];
669 /* copy values to argument arrays; set very small values to zero */
670 cpet[i]= va*cb[i] + vv*cvb + (1.0-vb)*(ct1+ct2);
671 if(fabs(cpet[i])<1.0e-12) cpet[i]=0.0;
672 if(cta!=NULL) {cta[i]=(1.0-vb)*ct1; if(fabs(cta[i])<1.0e-12) cta[i]=0.0;}
673 if(ctb!=NULL) {ctb[i]=(1.0-vb)*ct2; if(fabs(ctb[i])<1.0e-12) ctb[i]=0.0;}
674 if(ctab!=NULL) {ctab[i]=va*cb[i]; if(fabs(ctab[i])<1.0e-12) ctab[i]=0.0;}
675 if(ctvb!=NULL) {ctvb[i]=vv*cvb; if(fabs(ctvb[i])<1.0e-12) ctvb[i]=0.0;}
676 /* prepare to the next loop */
677 t_last=t[i]; ca_last=ca[i];
678 ct1_last=ct1; ct1i_last=ct1i;
679 ct2_last=ct2; ct2i_last=ct2i;
680 }
681
682 return 0;
683}
684/*****************************************************************************/
685
686/*****************************************************************************/
709 double *t,
711 double *ca,
713 double *cb,
715 int nr,
717 double k1,
719 double k2,
721 double k3,
723 double k4,
725 double k5,
727 double k6,
729 double kLoss,
731 double f,
733 double vb,
735 double fa,
737 double *cpet,
739 double *cta,
741 double *ctb,
743 double *ctc,
745 double *ctab,
747 double *ctvb
748) {
749 int i;
750 double dt2, b, c, d, w, z, u, va, vv;
751 double cai, ca_last, t_last, dct, cvb;
752 double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
753 double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
754
755
756 /* Check for data */
757 if(nr<2) return 1;
758 if(cpet==NULL) return 2;
759
760 /* Check parameters */
761 if(k1<0.0) return 3;
762 if(vb<0.0 || vb>=1.0) return 4;
763 if(fa<=0.0 || fa>1.0) return 5;
764 va=fa*vb; vv=(1.0-fa)*vb;
765
766 /* Calculate curves */
767 t_last=0.0; if(t[0]<t_last) t_last=t[0];
768 cai=ca_last=0.0;
769 ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=0.0;
770 ct1=ct2=ct3=ct1i=ct2i=ct3i=0.0;
771 for(i=0; i<nr; i++) {
772 /* delta time / 2 */
773 dt2=0.5*(t[i]-t_last);
774 /* calculate values */
775 if(dt2<0.0) {
776 return 5;
777 } else if(dt2>0.0) {
778 /* arterial integral */
779 cai+=(ca[i]+ca_last)*dt2;
780 /* Calculate partial results */
781 w = 1.0 + k4*dt2;
782 z = 1.0 + dt2*(k6 + kLoss);
783 u = k2 + k3 + k5 - k3*k4*dt2/w - k5*k6*dt2/z;
784 b = ct1i_last+dt2*ct1_last;
785 c = ct2i_last+dt2*ct2_last;
786 d = ct3i_last+dt2*ct3_last;
787 /* 1st tissue compartment and its integral */
788 ct1 = ( k1*cai - u*b + k4*c/w + k6*d/z ) / ( 1.0 + dt2*u);
789 ct1i = ct1i_last + dt2*(ct1_last+ct1);
790 /* 2nd tissue compartment and its integral */
791 ct2 = (k3*ct1i - k4*c) / w;
792 ct2i = ct2i_last + dt2*(ct2_last+ct2);
793 /* 3rd tissue compartment and its integral */
794 ct3 = (k5*ct1i - (k6 + kLoss)*d) / z;
795 ct3i = ct3i_last + dt2*(ct3_last+ct3);
796 }
797 /* Venous curve */
798 if(f>0.) {dct = k1*ca[i] - k2*ct1 - kLoss*ct3; cvb = cb[i] - dct/f;}
799 else cvb=cb[i];
800 /* copy values to argument arrays; set very small values to zero */
801 cpet[i]= va*cb[i] + vv*cvb + (1.0-vb)*(ct1+ct2+ct3);
802 if(fabs(cpet[i])<1.0e-12) cpet[i]=0.0;
803 if(cta!=NULL) {cta[i]=(1.0-vb)*ct1; if(fabs(cta[i])<1.0e-12) cta[i]=0.0;}
804 if(ctb!=NULL) {ctb[i]=(1.0-vb)*ct2; if(fabs(ctb[i])<1.0e-12) ctb[i]=0.0;}
805 if(ctc!=NULL) {ctc[i]=(1.0-vb)*ct3; if(fabs(ctc[i])<1.0e-12) ctc[i]=0.0;}
806 if(ctab!=NULL) {ctab[i]=va*cb[i]; if(fabs(ctab[i])<1.0e-12) ctab[i]=0.0;}
807 if(ctvb!=NULL) {ctvb[i]=vv*cvb; if(fabs(ctvb[i])<1.0e-12) ctvb[i]=0.0;}
808 /* prepare to the next loop */
809 t_last=t[i]; ca_last=ca[i];
810 ct1_last=ct1; ct1i_last=ct1i;
811 ct2_last=ct2; ct2i_last=ct2i;
812 ct3_last=ct3; ct3i_last=ct3i;
813 }
814
815 return 0;
816}
817/*****************************************************************************/
818
819/*****************************************************************************/
837 double *t,
839 double *cr,
841 int nr,
843 double R1,
845 double k2,
847 double k3,
849 double k4,
851 double *ct,
853 double *cta,
855 double *ctb
856) {
857 int i;
858 double f, b, w, dt2;
859 double cri, cr_last, t_last;
860 double cf, cf_last, cb, cb_last;
861 double cfi, cfi_last, cbi, cbi_last;
862
863
864 /* Check for data */
865 if(nr<2) return 1;
866 if(ct==NULL) return 2;
867
868 /* Calculate curves */
869 t_last=0.0; if(t[0]<t_last) t_last=t[0];
870 cri=cr_last=0.0; cf_last=cb_last=cfi_last=cbi_last=cf=cb=cfi=cbi=0.0;
871 for(i=0; i<nr; i++) {
872 /* delta time / 2 */
873 dt2=0.5*(t[i]-t_last);
874 /* calculate values */
875 if(dt2<0.0) {
876 return 5;
877 } else if(dt2>0.0) {
878 /* reference integral */
879 cri+=(cr[i]+cr_last)*dt2;
880 /* partial results */
881 f=cfi_last+dt2*cf_last;
882 b=cbi_last+dt2*cb_last;
883 w=k2 + k3 + k2*k4*dt2;
884 /* 1st tissue compartment and its integral */
885 cf = ( (1.0 + k4*dt2)*(R1*cr[i] + k2*cri) + k4*b - w*f ) /
886 ( 1.0 + dt2*(w+k4) );
887 cfi = cfi_last + dt2*(cf_last+cf);
888 /* 2nd tissue compartment and its integral */
889 cb = (k3*cfi - k4*b) / (1.0 + k4*dt2);
890 cbi = cbi_last + dt2*(cb_last+cb);
891 }
892 /* copy values to argument arrays; set very small values to zero */
893 ct[i]=cf+cb; if(fabs(ct[i])<1.0e-12) ct[i]=0.0;
894 if(cta!=NULL) {cta[i]=cf; if(fabs(cta[i])<1.0e-12) cta[i]=0.0;}
895 if(ctb!=NULL) {ctb[i]=cb; if(fabs(ctb[i])<1.0e-12) ctb[i]=0.0;}
896 /* prepare to the next loop */
897 t_last=t[i]; cr_last=cr[i];
898 cf_last=cf; cfi_last=cfi;
899 cb_last=cb; cbi_last=cbi;
900 }
901
902 return 0;
903}
904/*****************************************************************************/
905
906/*****************************************************************************/
921 double *t,
923 double *cr,
925 int nr,
927 double R1,
929 double k2,
931 double BP,
933 double *ct
934) {
935 int i;
936 double dt2;
937 double cri, cr_last, t_last;
938 double ct_last, cti, cti_last;
939
940
941 /* Check for data */
942 if(nr<2) return 1;
943 if(ct==NULL) return 2;
944
945 /* Calculate curves */
946 t_last=0.0; if(t[0]<t_last) t_last=t[0];
947 cri=cr_last=0.0; cti=ct_last=cti_last=0.0;
948 for(i=0; i<nr; i++) {
949 /* delta time / 2 */
950 dt2=0.5*(t[i]-t_last);
951 /* calculate values */
952 if(dt2<0.0) {
953 return 5;
954 } else if(dt2>0.0) {
955 /* reference integral */
956 cri+=(cr[i]+cr_last)*dt2;
957 /* Tissue compartment and its integral */
958 ct[i] = ( R1*cr[i] + k2*cri - (k2/(1.0+BP))*(cti_last+dt2*ct_last) ) /
959 ( 1.0 + dt2*(k2/(1.0+BP)) );
960 cti = cti_last + dt2*(ct_last+ct[i]);
961 }
962 /* set invalid or very small values to zero */
963 if(!(fabs(ct[i])>=1.0e-12)) ct[i]=0.0;
964 /* prepare to the next loop */
965 t_last=t[i]; cr_last=cr[i];
966 ct_last=ct[i]; cti_last=cti;
967 }
968
969 return 0;
970}
971/*****************************************************************************/
972
973/*****************************************************************************/
988 double *t,
990 double *cr,
992 int nr,
994 double R1,
996 double k2,
998 double k3,
1000 double *ct
1001) {
1002 int i;
1003 double dt2;
1004 double cri, cr_last, t_last;
1005 double ct_last, cti, cti_last;
1006
1007
1008 /* Check for data */
1009 if(nr<2) return 1;
1010 if(ct==NULL) return 2;
1011
1012 /* Calculate curves */
1013 t_last=0.0; if(t[0]<t_last) t_last=t[0];
1014 cri=cr_last=0.0; cti=ct_last=cti_last=0.0;
1015 for(i=0; i<nr; i++) {
1016 /* delta time / 2 */
1017 dt2=0.5*(t[i]-t_last);
1018 /* calculate values */
1019 if(dt2<0.0) {
1020 return 5;
1021 } else if(dt2>0.0) {
1022 /* reference integral */
1023 cri+=(cr[i]+cr_last)*dt2;
1024 /* Tissue compartment and its integral */
1025 ct[i] = ( R1*cr[i] + R1*k3*cri - (k2+k3)*(cti_last+dt2*ct_last) ) /
1026 ( 1.0 + dt2*(k2+k3) );
1027 cti = cti_last + dt2*(ct_last+ct[i]);
1028 }
1029 /* set very small values to zero */
1030 if(fabs(ct[i])<1.0e-12) ct[i]=0.0;
1031 /* prepare to the next loop */
1032 t_last=t[i]; cr_last=cr[i];
1033 ct_last=ct[i]; cti_last=cti;
1034 }
1035
1036 return 0;
1037}
1038/*****************************************************************************/
1039
1040/*****************************************************************************/
1056 double *t,
1058 double *ctot,
1060 int nr,
1062 double k01,
1064 double k12,
1066 double k21,
1068 double k03,
1070 double k34,
1072 double k43,
1074 double *c0,
1076 double *c1,
1078 double *c3
1079) {
1080 int i;
1081 double dt2, t_last, ictot, ctot_last;
1082 double c0_, c1_, c2_, c3_, c4_;
1083 double ic0_, ic1_, ic2_, ic3_, ic4_;
1084 double c0_last, c1_last, c2_last, c3_last, c4_last;
1085 double ic0_last, ic1_last, ic2_last, ic3_last, ic4_last;
1086 double a, b, am1, am2, am3, am4;
1087
1088
1089 /* Check input */
1090 if(t==NULL || ctot==NULL || nr<2) return(1);
1091 if(k01<0 || k12<0 || k21<0 || k03<0 || k34<0 || k43<0) return(2);
1092 if(t[0]<0.0) return(3);
1093
1094 /* Compute the TACs */
1095 t_last=0.0; if(t[0]<t_last) t_last=t[0];
1096 ictot=ctot_last=0.0;
1097 c0_=c0_last=ic0_=ic0_last=0.0;
1098 c1_=c1_last=ic1_=ic1_last=0.0;
1099 c2_=c2_last=ic2_=ic2_last=0.0;
1100 c3_=c3_last=ic3_=ic3_last=0.0;
1101 c4_=c4_last=ic4_=ic4_last=0.0;
1102 for(i=0; i<nr; i++) {
1103 /* delta time / 2 */
1104 dt2=0.5*(t[i]-t_last); if(dt2<0.0) return(5);
1105 if(dt2>0.0) {
1106 /* Compute temp constants */
1107 a=k01+k12-(k12*k21*dt2/(1.0+dt2*k21));
1108 b=k03+k34-(k34*k43*dt2/(1.0+dt2*k43));
1109 am1=ic1_last+dt2*c1_last;
1110 am2=ic2_last+dt2*c2_last;
1111 am3=ic3_last+dt2*c3_last;
1112 am4=ic4_last+dt2*c4_last;
1113
1114 /* Compute the "input" i.e. ctot integral */
1115 ictot+=(ctot[i]+ctot_last)*dt2;
1116 /* Compute C1(t) and its integral */
1117 c1_= ( k01*(1.0-k03*dt2/(1.0+dt2*b))*ictot
1118 -(a-k01*k03*dt2/(1.0+dt2*b))*am1
1119 +(k21/(1.0+dt2*k21))*am2
1120 -(k01/(1.0+dt2*b))*am3
1121 -(k01*k43*dt2/((1.0+dt2*b)*(1.0+dt2*k43)))*am4
1122 ) / ( 1.0+dt2*(a-k01*k03*dt2/(1.0+dt2*b)) );
1123 ic1_= ic1_last + dt2*(c1_+c1_last);
1124 /* Compute C2(t) and its integral */
1125 c2_= (k12*ic1_-k21*am2)/(1.0+dt2*k21);
1126 ic2_= ic2_last + dt2*(c2_+c2_last);
1127
1128 /* Compute C3(t) and its integral */
1129 c3_= ( k03*(1.0-k01*dt2/(1.0+dt2*a))*ictot
1130 -(b-k01*k03*dt2/(1.0+dt2*a))*am3
1131 +(k43/(1.0+dt2*k43))*am4
1132 -(k03/(1.0+dt2*a))*am1
1133 -(k03*k21*dt2/((1.0+dt2*a)*(1.0+dt2*k21)))*am2
1134 ) / ( 1.0+dt2*(b-k01*k03*dt2/(1.0+dt2*a)) );
1135 ic3_= ic3_last + dt2*(c3_+c3_last);
1136 /* Compute C4(t) and its integral */
1137 c4_= (k34*ic3_-k43*am4)/(1.0+dt2*k43);
1138 ic4_= ic4_last + dt2*(c4_+c4_last);
1139
1140 /* Compute the C0(t) */
1141 c0_=ctot[i]-c1_-c3_; /*if(c0_<0.0) return(-1);*/
1142 }
1143 /* Set output data */
1144 if(c0) c0[i]=c0_;
1145 if(c1) c1[i]=c1_;
1146 if(c3) c3[i]=c3_;
1147 /* Prepare for the next sample */
1148 c0_last=c0_; c1_last=c1_; c2_last=c2_; c3_last=c3_; c4_last=c4_;
1149 ic0_last=ic0_; ic1_last=ic1_; ic2_last=ic2_; ic3_last=ic3_; ic4_last=ic4_;
1150 ctot_last=ctot[i]; t_last=t[i];
1151 } /* next sample */
1152
1153 return(0);
1154}
1155/*****************************************************************************/
1156
1157/*****************************************************************************/
1167 double *t,
1169 double *ctot,
1171 int nr,
1173 double km,
1175 double k1m,
1177 double k2m,
1179 double k3m,
1181 double k4m,
1183 double *ca,
1185 double *cm
1186) {
1187 int i;
1188 double t_last, dt2;
1189 double ctoti, ctot_last;
1190 double a, b;
1191 double ct1m, ct1mi, ct1m_last, ct1mi_last;
1192 double ct2m, ct2mi, ct2m_last, ct2mi_last;
1193 double cpm, cpmi, cpm_last, cpmi_last;
1194
1195
1196 /* Check for data */
1197 if(nr<2) return 1;
1198 if(t==NULL || ctot==NULL || ca==NULL) return 2;
1199
1200 /* Calculate curves */
1201 t_last=0.0; if(t[0]<t_last) t_last=t[0];
1202 ctoti=ctot_last=0.0;
1203 ct1m_last=ct1mi_last=ct2m_last=ct2mi_last=cpm_last=cpmi_last=0.0;
1204 //ct1=ct2=ct3=ct1i=ct2i=ct3i=0.0;
1205 for(i=0; i<nr; i++) {
1206 /* delta time / 2 */
1207 dt2=0.5*(t[i]-t_last);
1208 /* calculate values */
1209 if(dt2<0.0) {
1210 return 5;
1211 } else {
1212 /* arterial integral */
1213 ctoti+=(ctot[i]+ctot_last)*dt2;
1214 /* partial results */
1215 a=k4m/(1.0+k4m*dt2);
1216 b=(k1m-km)/(1.0+dt2*(k1m+k3m-k3m*a*dt2));
1217 /* 1st tissue compartment and its integral */
1218 ct1m = (km*ctoti - k2m*(1.0-dt2*b)*(ct1mi_last+dt2*ct1m_last)
1219 + b*(cpmi_last+dt2*cpm_last) + a*b*dt2*(ct2mi_last+dt2*ct2m_last)
1220 ) / (1.0 + dt2*k2m*(1.0-dt2*b));
1221 ct1mi = ct1mi_last + dt2*(ct1m_last+ct1m);
1222 /* Metabolite plasma compartment and its integral */
1223 cpm = (k2m*ct1mi + a*(ct2mi_last+dt2*ct2m_last)
1224 - (k1m+k3m-k3m*dt2*a)*(cpmi_last+dt2*cpm_last)
1225 ) / (1.0 + dt2*(k1m+k3m-k3m*dt2*a));
1226 cpmi = cpmi_last + dt2*(cpm_last+cpm);
1227 /* 2nd tissue compartment and its integral */
1228 ct2m = (k3m*cpmi - k4m*(ct2mi_last+dt2*ct2m_last)) / (1.0 + dt2*k4m);
1229 ct2mi = ct2mi_last + dt2*(ct2m_last+ct2m);
1230 }
1231 /* copy values to argument arrays; set very small values to zero */
1232 if(ca!=NULL) {ca[i]=ctot[i]-cpm; if(fabs(ca[i])<1.0e-12) ca[i]=0.0;}
1233 if(cm!=NULL) {cm[i]=cpm; if(fabs(cm[i])<1.0e-12) cm[i]=0.0;}
1234 /* prepare to the next loop */
1235 t_last=t[i]; ctot_last=ctot[i];
1236 ct1m_last=ct1m; ct1mi_last=ct1mi;
1237 ct2m_last=ct2m; ct2mi_last=ct2mi;
1238 cpm_last=cpm; cpmi_last=cpmi;
1239 } // next sample
1240 return 0;
1241}
1242/*****************************************************************************/
1243
1244/*****************************************************************************/
1255 double *t,
1257 double *ci,
1259 int nr,
1261 double k1,
1263 double k2,
1265 double Vfit,
1267 double *ct
1268) {
1269 int i;
1270 double dt2;
1271 double cii, ci_last, t_last;
1272 double ct_last, cti, cti_last;
1273
1274
1275 /* Check for data */
1276 if(nr<2) return(1);
1277 if(ct==NULL) return(2);
1278
1279 /* Calculate curves */
1280 t_last=0.0; cii=ci_last=0.0;
1281 cti=ct_last=cti_last=0.0;
1282 for(i=0; i<nr; i++) {
1283 /* delta time / 2 */
1284 dt2=0.5*(t[i]-t_last);
1285 /* calculate values */
1286 if(dt2<0.0) {
1287 return(5);
1288 } else if(dt2>0.0) {
1289 /* input integral */
1290 cii+=(ci[i]+ci_last)*dt2;
1291 /* Tissue compartment and its integral */
1292 ct[i] = (Vfit*ci[i] + k1*cii - k2*(cti_last+dt2*ct_last)) / (1.0 + dt2*k2);
1293 cti = cti_last + dt2*(ct_last+ct[i]);
1294 }
1295 /* set very small values to zero */
1296 if(fabs(ct[i])<1.0e-12) ct[i]=0.0;
1297 /* prepare to the next loop */
1298 t_last=t[i]; ci_last=ci[i];
1299 ct_last=ct[i]; cti_last=cti;
1300 }
1301 return(0);
1302}
1303/*****************************************************************************/
1304
1305/*****************************************************************************/
1319 double *t,
1321 double *ca,
1323 int nr,
1325 double k1,
1327 double k2,
1329 double *ct
1330) {
1331 int i;
1332 double dt2;
1333 double cai, ca_last, t_last;
1334 double ct1, ct1_last;
1335 double ct1i, ct1i_last;
1336
1337
1338 /* Check for data */
1339 if(nr<2) return 1;
1340 if(ct==NULL) return 2;
1341
1342 /* Check actual parameter number */
1343 if(k1<0.0) return 3;
1344
1345 /* Calculate curves */
1346 t_last=0.0; if(t[0]<t_last) t_last=t[0];
1347 cai=ca_last=0.0;
1348 ct1_last=ct1i_last=0.0;
1349 ct1=ct1i=0.0;
1350 for(i=0; i<nr; i++) {
1351 /* delta time / 2 */
1352 dt2=0.5*(t[i]-t_last);
1353 /* calculate values */
1354 if(dt2<0.0) {
1355 return 5;
1356 } else if(dt2>0.0) {
1357 /* arterial integral */
1358 cai+=(ca[i]+ca_last)*dt2;
1359 /* tissue compartment and its integral */
1360 ct1 = (k1*cai - k2*(ct1i_last+dt2*ct1_last)) / (1.0 + dt2*k2);
1361 ct1i = ct1i_last + dt2*(ct1_last+ct1);
1362 }
1363 /* copy values to argument arrays; set very small values to zero */
1364 ct[i]=ct1; if(fabs(ct[i])<1.0e-12) ct[i]=0.0;
1365 /* prepare to the next loop */
1366 t_last=t[i]; ca_last=ca[i];
1367 ct1_last=ct1; ct1i_last=ct1i;
1368 }
1369
1370 return 0;
1371}
1372/*****************************************************************************/
1373
1374/*****************************************************************************/
1389 double *t,
1391 double *ca1,
1393 double *ca2,
1395 double *cb,
1397 int nr,
1399 double k1,
1401 double k2,
1403 double k3,
1405 double k4,
1407 double k5,
1409 double k6,
1411 double k1b,
1413 double k2b,
1415 double f,
1417 double vb,
1419 double fa,
1421 double *scpet,
1423 double *sct1,
1425 double *sct2,
1427 double *sct3,
1429 double *sct1b,
1431 double *sctab,
1433 double *sctvb
1434) {
1435 int i;
1436 double b, c, d, e, w, z, dt2, va, vv;
1437 double ca1i, ca1_last, ca2i, ca2_last, t_last, dct, cvb;
1438 double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
1439 double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
1440 double ct1b, ct1b_last, ct1bi, ct1bi_last;
1441
1442
1443 /* Check for data */
1444 if(nr<2) return 1;
1445 if(scpet==NULL) return 2;
1446
1447 /* Check parameters */
1448 if(k1<0.0) return 3;
1449 if(vb<0.0 || vb>=1.0) return 4;
1450 if(fa<=0.0 || fa>1.0) return 5;
1451 va=fa*vb; vv=(1.0-fa)*vb;
1452
1453 /* Calculate curves */
1454 t_last=0.0; if(t[0]<t_last) t_last=t[0];
1455 ca1i=ca1_last=ca2i=ca2_last=0.0;
1456 ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=0.0;
1457 ct1b_last=ct1bi_last=0.0;
1458 ct1=ct2=ct3=ct1b=ct1i=ct2i=ct3i=ct1bi=0.0;
1459 for(i=0; i<nr; i++) {
1460 /* delta time / 2 */
1461 dt2=0.5*(t[i]-t_last);
1462 /* calculate values */
1463 if(dt2<0.0) {
1464 return 5;
1465 } else if(dt2>0.0) {
1466 /* arterial integrals */
1467 ca1i+=(ca1[i]+ca1_last)*dt2;
1468 ca2i+=(ca2[i]+ca2_last)*dt2;
1469 /* partial results */
1470 b=ct1i_last+dt2*ct1_last;
1471 c=ct2i_last+dt2*ct2_last;
1472 d=ct3i_last+dt2*ct3_last;
1473 e=ct1bi_last+dt2*ct1b_last;
1474 w=k4 + k5 - (k5*k6*dt2)/(1.0+k6*dt2);
1475 z=1.0+w*dt2;
1476 /* 1st tissue compartment and its integral */
1477 ct1 = (
1478 + k1*z*ca1i + (k3*k4*dt2 - (k2+k3)*z)*b
1479 + k4*c + k4*k6*dt2*d/(1.0+k6*dt2)
1480 ) / ( z*(1.0 + dt2*(k2+k3)) - k3*k4*dt2*dt2 );
1481 ct1i = ct1i_last + dt2*(ct1_last+ct1);
1482 ct1b = (k1b*ca2i - k2b*e) / (1.0 + dt2*k2b);
1483 ct1bi = ct1bi_last + dt2*(ct1b_last+ct1b);
1484 /* 2nd tissue compartment and its integral */
1485 ct2 = (k3*ct1i - w*c + k6*d/(1.0+k6*dt2)) / z;
1486 ct2i = ct2i_last + dt2*(ct2_last+ct2);
1487 /* 3rd tissue compartment and its integral */
1488 ct3 = (k5*ct2i - k6*d) / (1.0 + k6*dt2);
1489 ct3i = ct3i_last + dt2*(ct3_last+ct3);
1490 }
1491 /* Venous curve */
1492 if(f>0.) {
1493 dct = k1*ca1[i] - k2*ct1 + k1b*ca2[i] - k2b*ct1b;
1494 cvb = cb[i] - dct/f;
1495 } else cvb=cb[i];
1496 /* copy values to argument arrays; set very small values to zero */
1497 scpet[i]= va*cb[i] + vv*cvb + (1.0-vb)*(ct1+ct2+ct3+ct1b);
1498 if(fabs(scpet[i])<1.0e-12) scpet[i]=0.0;
1499 if(sct1!=NULL) {sct1[i]=(1.0-vb)*ct1; if(fabs(sct1[i])<1.0e-12) sct1[i]=0.0;}
1500 if(sct2!=NULL) {sct2[i]=(1.0-vb)*ct2; if(fabs(sct2[i])<1.0e-12) sct2[i]=0.0;}
1501 if(sct3!=NULL) {sct3[i]=(1.0-vb)*ct3; if(fabs(sct3[i])<1.0e-12) sct3[i]=0.0;}
1502 if(sct1b!=NULL) {
1503 sct1b[i]=(1.0-vb)*ct1b; if(fabs(sct1b[i])<1.0e-12) sct1b[i]=0.0;}
1504 if(sctab!=NULL) {sctab[i]=va*cb[i]; if(fabs(sctab[i])<1.0e-12) sctab[i]=0.0;}
1505 if(sctvb!=NULL) {sctvb[i]=vv*cvb; if(fabs(sctvb[i])<1.0e-12) sctvb[i]=0.0;}
1506 /* prepare to the next loop */
1507 t_last=t[i]; ca1_last=ca1[i]; ca2_last=ca2[i];
1508 ct1_last=ct1; ct1i_last=ct1i;
1509 ct2_last=ct2; ct2i_last=ct2i;
1510 ct3_last=ct3; ct3i_last=ct3i;
1511 ct1b_last=ct1b; ct1bi_last=ct1bi;
1512 }
1513
1514 return 0;
1515}
1516/*****************************************************************************/
1517
1518/*****************************************************************************/
1535 double *t,
1537 double *ca1,
1539 double *ca2,
1541 double *cb,
1543 int nr,
1545 double k1,
1547 double k2,
1549 double k3,
1551 double k4,
1553 double k5,
1555 double k6,
1557 double k7,
1559 double km,
1561 double k1b,
1563 double k2b,
1565 double f,
1567 double vb,
1569 double fa,
1571 double *scpet,
1573 double *sct1,
1575 double *sct2,
1577 double *sct3,
1579 double *sct1b,
1581 double *sctab,
1583 double *sctvb,
1585 int verbose
1586) {
1587 int i;
1588 double b, c, d, e, pt, qt, dt2, va, vv;
1589 double ca1i, ca1_last, ca2i, ca2_last, t_last, dct, cvb;
1590 double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
1591 double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
1592 double ct1b, ct1b_last, ct1bi, ct1bi_last;
1593
1594
1595 if(verbose>0) {
1596 printf("simC4DIvp()\n");
1597 if(verbose>1) {
1598 printf(" k1 := %g\n", k1);
1599 printf(" k2 := %g\n", k2);
1600 printf(" k3 := %g\n", k3);
1601 printf(" k4 := %g\n", k4);
1602 printf(" k5 := %g\n", k5);
1603 printf(" k6 := %g\n", k6);
1604 printf(" k7 := %g\n", k7);
1605 printf(" km := %g\n", km);
1606 printf(" k1b := %g\n", k1b);
1607 printf(" k2b := %g\n", k2b);
1608 printf(" vb := %g\n", vb);
1609 printf(" fa := %g\n", fa);
1610 printf(" f := %g\n", f);
1611 }
1612 }
1613
1614 /* Check for data */
1615 if(nr<2) return 1;
1616 if(scpet==NULL) return 2;
1617
1618 /* Check parameters */
1619 if(k1<0.0 || k1b<0.0) return 3;
1620 if(vb<0.0 || vb>=1.0) return 4;
1621 if(fa<0.0 || fa>1.0) return 5;
1622 va=fa*vb; vv=(1.0-fa)*vb;
1623
1624 /* Calculate curves */
1625 t_last=0.0; if(t[0]<t_last) t_last=t[0];
1626 ca1i=ca1_last=ca2i=ca2_last=0.0;
1627 ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=ct1b_last=
1628 ct1bi_last=0.0;
1629 ct1=ct2=ct3=ct1b=ct1i=ct2i=ct3i=ct1bi=0.0;
1630 for(i=0; i<nr; i++) {
1631 /* delta time / 2 */
1632 dt2=0.5*(t[i]-t_last);
1633 /* calculate values */
1634 if(dt2<0.0) {
1635 return 5;
1636 } else if(dt2>0.0) {
1637 /* arterial integrals */
1638 ca1i+=(ca1[i]+ca1_last)*dt2;
1639 ca2i+=(ca2[i]+ca2_last)*dt2;
1640 /* partial results */
1641 b=ct1i_last+dt2*ct1_last;
1642 c=ct2i_last+dt2*ct2_last;
1643 d=ct3i_last+dt2*ct3_last;
1644 e=ct1bi_last+dt2*ct1b_last;
1645 pt=k6+k7;
1646 qt=k2+k3+k5+km-(k3*k4*dt2)/(1.0+k4*dt2)-(k5*k6*dt2)/(1.0+pt*dt2);
1647 /* 1st tissue compartment and its integral */
1648 ct1 = (k1/(1.0+qt*dt2))*ca1i
1649 - (qt/(1.0+qt*dt2))*b
1650 + (k4/((1.0+qt*dt2)*(1.0+k4*dt2)))*c
1651 + (k6/((1.0+qt*dt2)*(1.0+pt*dt2)))*d;
1652 ct1i = ct1i_last + dt2*(ct1_last+ct1);
1653 /* 2nd tissue compartment and its integral */
1654 ct2 = (k3/(1.0+k4*dt2))*ct1i
1655 - (k4/(1.0+k4*dt2))*c;
1656 ct2i = ct2i_last + dt2*(ct2_last+ct2);
1657 /* 3rd tissue compartment and its integral */
1658 ct3 = (k5/(1.0+pt*dt2))*ct1i
1659 - (pt/(1.0+pt*dt2))*d;
1660 ct3i = ct3i_last + dt2*(ct3_last+ct3);
1661 /* 4th tissue compartment (the 1st for tracer 2) and its integral */
1662 ct1b = (k1b/(1.0+k2b*dt2))*ca2i
1663 - (k2b/(1.0+k2b*dt2))*e
1664 + (km/(1.0+k2b*dt2))*ct1i;
1665 ct1bi = ct1bi_last + dt2*(ct1b_last+ct1b);
1666 }
1667 /* Venous curve */
1668 if(f>0.) {
1669 dct = k1*ca1[i] - k2*ct1 - k7*ct3 + k1b*ca2[i] - k2b*ct1b;
1670 cvb = cb[i] - dct/f;
1671 } else cvb=cb[i];
1672 /* copy values to argument arrays; set very small values to zero */
1673 scpet[i]= va*cb[i] + vv*cvb + (1.0-vb)*(ct1+ct2+ct3+ct1b);
1674 if(fabs(scpet[i])<1.0e-12) scpet[i]=0.0;
1675 if(sct1!=NULL) {
1676 sct1[i]=(1.0-vb)*ct1; if(fabs(sct1[i])<1.0e-12) sct1[i]=0.0;}
1677 if(sct2!=NULL) {
1678 sct2[i]=(1.0-vb)*ct2; if(fabs(sct2[i])<1.0e-12) sct2[i]=0.0;}
1679 if(sct3!=NULL) {
1680 sct3[i]=(1.0-vb)*ct3; if(fabs(sct3[i])<1.0e-12) sct3[i]=0.0;}
1681 if(sct1b!=NULL) {
1682 sct1b[i]=(1.0-vb)*ct1b; if(fabs(sct1b[i])<1.0e-12) sct1b[i]=0.0;}
1683 if(sctab!=NULL) {
1684 sctab[i]=va*cb[i]; if(fabs(sctab[i])<1.0e-12) sctab[i]=0.0;}
1685 if(sctvb!=NULL) {
1686 sctvb[i]=vv*cvb; if(fabs(sctvb[i])<1.0e-12) sctvb[i]=0.0;}
1687 /* prepare to the next loop */
1688 t_last=t[i]; ca1_last=ca1[i]; ca2_last=ca2[i];
1689 ct1_last=ct1; ct1i_last=ct1i;
1690 ct2_last=ct2; ct2i_last=ct2i;
1691 ct3_last=ct3; ct3i_last=ct3i;
1692 ct1b_last=ct1b; ct1bi_last=ct1bi;
1693 }
1694
1695 if(verbose>2) {
1696 printf("AUC 0-%g:\n", t_last);
1697 printf(" ca1i := %g\n", ca1i);
1698 printf(" ca2i := %g\n", ca2i);
1699 printf(" ct1i := %g\n", ct1i_last);
1700 printf(" ct2i := %g\n", ct2i_last);
1701 printf(" ct3i := %g\n", ct3i_last);
1702 printf(" ct1bi := %g\n", ct1bi_last);
1703 }
1704
1705 return 0;
1706}
1707/*****************************************************************************/
1708
1709/*****************************************************************************/
1726 double *t,
1728 double *ca1,
1730 double *ca2,
1732 double *cb,
1734 int nr,
1736 double k1,
1738 double k2,
1740 double k3,
1742 double k4,
1744 double k5,
1746 double k6,
1748 double k7,
1750 double km,
1752 double k1b,
1754 double k2b,
1756 double f,
1758 double vb,
1760 double fa,
1762 double *scpet,
1764 double *sct1,
1766 double *sct2,
1768 double *sct3,
1770 double *sct1b,
1772 double *sctab,
1774 double *sctvb,
1776 int verbose
1777) {
1778 int i;
1779 double b, c, d, e, pt, qt, rt, dt2, va, vv;
1780 double ca1i, ca1_last, ca2i, ca2_last, t_last, dct, cvb;
1781 double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
1782 double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
1783 double ct1b, ct1b_last, ct1bi, ct1bi_last;
1784
1785
1786 if(verbose>0) {
1787 printf("simC4DIvs()\n");
1788 if(verbose>1) {
1789 printf(" k1 := %g\n", k1);
1790 printf(" k2 := %g\n", k2);
1791 printf(" k3 := %g\n", k3);
1792 printf(" k4 := %g\n", k4);
1793 printf(" k5 := %g\n", k5);
1794 printf(" k6 := %g\n", k6);
1795 printf(" k7 := %g\n", k7);
1796 printf(" km := %g\n", km);
1797 printf(" k1b := %g\n", k1b);
1798 printf(" k2b := %g\n", k2b);
1799 printf(" vb := %g\n", vb);
1800 printf(" fa := %g\n", fa);
1801 printf(" f := %g\n", f);
1802 }
1803 }
1804
1805 /* Check for data */
1806 if(nr<2) return 1;
1807 if(scpet==NULL) return 2;
1808
1809 /* Check parameters */
1810 if(k1<0.0 || k1b<0.0) return 3;
1811 if(vb<0.0 || vb>=1.0) return 4;
1812 if(fa<0.0 || fa>1.0) return 5;
1813 va=fa*vb; vv=(1.0-fa)*vb;
1814
1815 /* Calculate curves */
1816 t_last=0.0; if(t[0]<t_last) t_last=t[0];
1817 ca1i=ca1_last=ca2i=ca2_last=0.0;
1818 ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=ct1b_last=
1819 ct1bi_last=0.0;
1820 ct1=ct2=ct3=ct1b=ct1i=ct2i=ct3i=ct1bi=0.0;
1821 for(i=0; i<nr; i++) {
1822 /* delta time / 2 */
1823 dt2=0.5*(t[i]-t_last);
1824 /* calculate values */
1825 if(dt2<0.0) {
1826 return 5;
1827 } else if(dt2>0.0) {
1828 /* arterial integrals */
1829 ca1i+=(ca1[i]+ca1_last)*dt2;
1830 ca2i+=(ca2[i]+ca2_last)*dt2;
1831 //printf("ca1[%d]=%g int=%g\n", i, ca1[i], ca1i);
1832 //printf("ca2[%d]=%g int=%g\n", i, ca2[i], ca2i);
1833 /* partial results */
1834 b=ct1i_last+dt2*ct1_last;
1835 c=ct2i_last+dt2*ct2_last;
1836 d=ct3i_last+dt2*ct3_last;
1837 e=ct1bi_last+dt2*ct1b_last;
1838 pt=k6+k7;
1839 qt=k4+k5-(k5*k6*dt2)/(1.0+pt*dt2);
1840 rt=k2+k3+km-(k3*k4*dt2)/(1.0+qt*dt2);
1841 /* 1st tissue compartment and its integral */
1842 ct1 = (k1/(1.0+rt*dt2))*ca1i
1843 - (rt/(1.0+rt*dt2))*b
1844 + (k4/((1.0+qt*dt2)*(1.0+rt*dt2)))*c
1845 + ((k4*k6*dt2)/((1.0+pt*dt2)*(1.0+qt*dt2)*(1.0+rt*dt2)))*d;
1846 ct1i = ct1i_last + dt2*(ct1_last+ct1);
1847 /* 2nd tissue compartment and its integral */
1848 ct2 = (k3/(1.0+qt*dt2))*ct1i
1849 - (qt/(1.0+qt*dt2))*c
1850 + (k6/((1.0+pt*dt2)*(1.0+qt*dt2)))*d;
1851 ct2i = ct2i_last + dt2*(ct2_last+ct2);
1852 /* 3rd tissue compartment and its integral */
1853 ct3 = (k5/(1.0+pt*dt2))*ct2i
1854 - (pt/(1.0+pt*dt2))*d;
1855 ct3i = ct3i_last + dt2*(ct3_last+ct3);
1856 /* 4th tissue compartment (the 1st for tracer 2) and its integral */
1857 ct1b = (k1b/(1.0+k2b*dt2))*ca2i
1858 - (k2b/(1.0+k2b*dt2))*e
1859 + (km/(1.0+k2b*dt2))*ct1i;
1860 ct1bi = ct1bi_last + dt2*(ct1b_last+ct1b);
1861 }
1862 /* Venous curve */
1863 if(f>0.) {
1864 dct = k1*ca1[i] - k2*ct1 - k7*ct3 + k1b*ca2[i] - k2b*ct1b;
1865 cvb = cb[i] - dct/f;
1866 } else cvb=cb[i];
1867 /* copy values to argument arrays; set very small values to zero */
1868 scpet[i]= va*cb[i] + vv*cvb + (1.0-vb)*(ct1+ct2+ct3+ct1b);
1869 if(fabs(scpet[i])<1.0e-12) scpet[i]=0.0;
1870 if(sct1!=NULL) {
1871 sct1[i]=(1.0-vb)*ct1; if(fabs(sct1[i])<1.0e-12) sct1[i]=0.0;}
1872 if(sct2!=NULL) {
1873 sct2[i]=(1.0-vb)*ct2; if(fabs(sct2[i])<1.0e-12) sct2[i]=0.0;}
1874 if(sct3!=NULL) {
1875 sct3[i]=(1.0-vb)*ct3; if(fabs(sct3[i])<1.0e-12) sct3[i]=0.0;}
1876 if(sct1b!=NULL) {
1877 sct1b[i]=(1.0-vb)*ct1b; if(fabs(sct1b[i])<1.0e-12) sct1b[i]=0.0;}
1878 if(sctab!=NULL) {
1879 sctab[i]=va*cb[i]; if(fabs(sctab[i])<1.0e-12) sctab[i]=0.0;}
1880 if(sctvb!=NULL) {
1881 sctvb[i]=vv*cvb; if(fabs(sctvb[i])<1.0e-12) sctvb[i]=0.0;}
1882 /* prepare to the next loop */
1883 t_last=t[i]; ca1_last=ca1[i]; ca2_last=ca2[i];
1884 ct1_last=ct1; ct1i_last=ct1i;
1885 ct2_last=ct2; ct2i_last=ct2i;
1886 ct3_last=ct3; ct3i_last=ct3i;
1887 ct1b_last=ct1b; ct1bi_last=ct1bi;
1888 }
1889
1890 if(verbose>2) {
1891 printf("AUC 0-%g:\n", t_last);
1892 printf(" ca1i := %g\n", ca1i);
1893 printf(" ca2i := %g\n", ca2i);
1894 printf(" ct1i := %g\n", ct1i_last);
1895 printf(" ct2i := %g\n", ct2i_last);
1896 printf(" ct3i := %g\n", ct3i_last);
1897 printf(" ct1bi := %g\n", ct1bi_last);
1898 }
1899
1900 return 0;
1901}
1902/*****************************************************************************/
1903
1904/*****************************************************************************/
1913 double *x,
1915 double *y,
1917 int n,
1919 double tau1,
1921 double tau2,
1924 double *tmp
1925) {
1926 /* Check input */
1927 if(x==NULL || y==NULL || n<2) return 1;
1928 if(tau1<0.0 || tau2<0.0) return 2;
1929
1930 /* Allocate memory if not allocated by user */
1931 double *buf;
1932 if(tmp!=NULL) buf=tmp; else buf=(double*)malloc(n*sizeof(double));
1933 if(buf==NULL) return 3;
1934
1935 /* First dispersion */
1936 if(tau1>0.0) {
1937 double k=1.0/tau1;
1938 int ret=simC1(x, y, n, k, k, buf);
1939 if(ret!=0) {
1940 if(tmp==NULL) free(buf);
1941 return 100+ret;
1942 }
1943 for(int i=0; i<n; i++) y[i]=buf[i];
1944 }
1945
1946 /* Second dispersion */
1947 if(tau2>0.0) {
1948 double k=1.0/tau2;
1949 int ret=simC1(x, y, n, k, k, buf);
1950 if(ret!=0) {
1951 if(tmp==NULL) free(buf);
1952 return 200+ret;
1953 }
1954 for(int i=0; i<n; i++) y[i]=buf[i];
1955 }
1956
1957 if(tmp==NULL) free(buf);
1958 return 0;
1959}
1960/*****************************************************************************/
1961
1962/*****************************************************************************/
1980 double *t,
1982 double *ca1,
1984 double *ca2,
1986 double *ca1i,
1988 double *ca2i,
1990 const int n,
1992 const double k1a,
1994 const double k2a,
1996 const double km,
1998 const double k1b,
2000 const double k2b,
2002 const double vb,
2004 const double fa,
2007 double *scpet,
2010 double *sct1,
2013 double *sct2,
2016 double *sctab,
2019 double *sctvb1,
2022 double *sctvb2,
2024 double *scvb1,
2026 double *scvb2,
2028 const int verbose
2029) {
2030 if(verbose>0) {
2031 printf("simOxygen()\n");
2032 if(verbose>1) {
2033 printf(" k1a := %g\n", k1a);
2034 printf(" k2a := %g\n", k2a);
2035 printf(" km := %g\n", km);
2036 printf(" k1b := %g\n", k1b);
2037 printf(" k2b := %g\n", k2b);
2038 printf(" vb := %g\n", vb);
2039 printf(" fa := %g\n", fa);
2040 printf(" n := %d\n", n);
2041 }
2042 }
2043
2044 /* Check for data */
2045 if(n<2) return 1;
2046 if(t==NULL) return 2;
2047 if(ca1==NULL || ca2==NULL) return 3;
2048
2049 /* Check parameters */
2050 if(k1a<0.0 || k1b<0.0 || k2a<0.0 || k2b<0.0) return(4);
2051 if(vb<0.0 || vb>=1.0) return(5);
2052 if(fa<0.0 || fa>1.0) return(6);
2053 double va=fa*vb; // arterial volume fraction in tissue
2054 double vv=(1.0-fa)*vb; // venous volume fraction in tissue
2055
2056 /* Set initial condition */
2057 double t_last=0.0; // previous sample time
2058 if(t[0]<t_last) t_last=t[0];
2059 /* Concentrations, integrals, and previous values are zero */
2060 double cba1i=0.0, cba1_last=0.0;
2061 double cba2i=0.0, cba2_last=0.0;
2062 double ct1=0.0, ct1_last=0.0, ct1i=0.0, ct1i_last=0.0;
2063 double ct2=0.0, ct2_last=0.0, ct2i=0.0, ct2i_last=0.0;
2064 double cvb1=0.0, cvb2=0.0;
2065
2066 /* Calculate curves */
2067 double p, q;
2068 double dt2; // delta t / 2
2069 for(int i=0; i<n; i++) {
2070 dt2=0.5*(t[i]-t_last);
2071 if(dt2>0.0) {
2072 /* arterial integrals */
2073 if(ca1i!=NULL) cba1i=ca1i[i]; else cba1i+=(ca1[i]+cba1_last)*dt2;
2074 if(ca2i!=NULL) cba2i=ca2i[i]; else cba2i+=(ca2[i]+cba2_last)*dt2;
2075 /* partial results */
2076 p=ct1i_last+dt2*ct1_last;
2077 q=ct2i_last+dt2*ct2_last;
2078 /* 1st tissue compartment and its integral */
2079 ct1 = (k1a*cba1i - (k2a+km)*p) / (1.0 + dt2*(k2a+km));
2080 ct1i = ct1i_last + dt2*(ct1_last+ct1);
2081 /* 2nd tissue compartment and its integral */
2082 ct2 = (km*ct1i + k1b*cba2i - k2b*q) / (1.0 + dt2*k2b);
2083 ct2i = ct2i_last + dt2*(ct2_last+ct2);
2084 }
2085 /* Venous BTACs */
2086 if(k1a>0.0 && k2a>0.0) cvb1=ct1/(k1a/k2a);
2087 else if(k2a>0.0) cvb1=0.0; else cvb1=ca1[i];
2088 if(k1b>0.0 && k2b>0.0) cvb2=ct2/(k1b/k2b);
2089 else if(k2b>0.0) cvb2=0.0; else cvb2=ca2[i];
2090 /* copy values to argument arrays; set very small values to zero */
2091 if(scpet!=NULL) {
2092 scpet[i]= va*(ca1[i]+ca2[i]) + vv*(cvb1+cvb2) + (1.0-vb)*(ct1+ct2);
2093 if(fabs(scpet[i])<1.0e-12) scpet[i]=0.0;
2094 }
2095 if(sct1!=NULL) {
2096 sct1[i]=(1.0-vb)*ct1;
2097 if(fabs(sct1[i])<1.0e-12) sct1[i]=0.0;
2098 }
2099 if(sct2!=NULL) {
2100 sct2[i]=(1.0-vb)*ct2;
2101 if(fabs(sct2[i])<1.0e-12) sct2[i]=0.0;
2102 }
2103 if(sctab!=NULL) {
2104 sctab[i]=va*(ca1[i]+ca2[i]);
2105 if(fabs(sctab[i])<1.0e-12) sctab[i]=0.0;
2106 }
2107 if(sctvb1!=NULL) {
2108 sctvb1[i]=vv*cvb1;
2109 if(fabs(sctvb1[i])<1.0e-12) sctvb1[i]=0.0;
2110 }
2111 if(sctvb2!=NULL) {
2112 sctvb2[i]=vv*cvb2;
2113 if(fabs(sctvb2[i])<1.0e-12) sctvb2[i]=0.0;
2114 }
2115 if(scvb1!=NULL) scvb1[i]=cvb1;
2116 if(scvb2!=NULL) scvb2[i]=cvb2;
2117 /* prepare for the next loop */
2118 t_last=t[i];
2119 cba1_last=ca1[i]; cba2_last=ca2[i];
2120 ct1_last=ct1; ct1i_last=ct1i;
2121 ct2_last=ct2; ct2i_last=ct2i;
2122 } // next sample
2123
2124 if(verbose>2) {
2125 printf("AUC 0-%g:\n", t_last);
2126 printf(" cba1i := %g\n", cba1i);
2127 printf(" cba2i := %g\n", cba2i);
2128 printf(" ct1i := %g\n", ct1i_last);
2129 printf(" ct2i := %g\n", ct2i_last);
2130 }
2131
2132 return 0;
2133}
2134/*****************************************************************************/
2135
2136/*****************************************************************************/
Header file for libtpcmodel.
int simTRTM(double *t, double *cr, int nr, double R1, double k2, double k3, double *ct)
Definition simulate.c:986
int simC3vp(double *t, double *ca, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double f, double vb, double fa, double *cpet, double *cta, double *ctb, double *ctc, double *ctab, double *ctvb)
Definition simulate.c:373
int simC4DIvp(double *t, double *ca1, double *ca2, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double k7, double km, double k1b, double k2b, double f, double vb, double fa, double *scpet, double *sct1, double *sct2, double *sct3, double *sct1b, double *sctab, double *sctvb, int verbose)
Definition simulate.c:1533
int simTPCMOD0009c(double *t, double *ctot, int nr, double km, double k1m, double k2m, double k3m, double k4m, double *ca, double *cm)
Definition simulate.c:1165
int simOxygen(double *t, double *ca1, double *ca2, double *ca1i, double *ca2i, const int n, const double k1a, const double k2a, const double km, const double k1b, const double k2b, const double vb, const double fa, double *scpet, double *sct1, double *sct2, double *sctab, double *sctvb1, double *sctvb2, double *scvb1, double *scvb2, const int verbose)
Definition simulate.c:1978
int simHuangmet(double *t, double *ctot, int nr, double k01, double k12, double k21, double k03, double k34, double k43, double *c0, double *c1, double *c3)
Definition simulate.c:1054
int simC3vs(double *t, double *ca, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double f, double vb, double fa, double *cpet, double *cta, double *ctb, double *ctc, double *ctab, double *ctvb)
Definition simulate.c:243
int simDispersion(double *x, double *y, int n, double tau1, double tau2, double *tmp)
Definition simulate.c:1911
int simC4DIvs(double *t, double *ca1, double *ca2, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double k7, double km, double k1b, double k2b, double f, double vb, double fa, double *scpet, double *sct1, double *sct2, double *sct3, double *sct1b, double *sctab, double *sctvb, int verbose)
Definition simulate.c:1724
int simSRTM(double *t, double *cr, int nr, double R1, double k2, double BP, double *ct)
Definition simulate.c:919
int simC1(double *t, double *ca, int nr, double k1, double k2, double *ct)
Definition simulate.c:1317
int simRTCM(double *t, double *cr, int nr, double R1, double k2, double k3, double k4, double *ct, double *cta, double *ctb)
Definition simulate.c:835
int simMBF(double *t, double *ci, int nr, double k1, double k2, double Vfit, double *ct)
Definition simulate.c:1253
int simC3vpKLoss(double *t, double *ca, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double kLoss, double f, double vb, double fa, double *cpet, double *cta, double *ctb, double *ctc, double *ctab, double *ctvb)
Definition simulate.c:707
int simC2vl(double *t, double *ca, double *cb, int nr, double k1, double k2, double k3, double kL, double f, double vb, double fa, double *cpet, double *cta, double *ctb, double *ctab, double *ctvb)
Definition simulate.c:592
int simC3p(double *t, double *ca, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double *ct, double *cta, double *ctb, double *ctc)
Definition simulate.c:136
int simC2l(double *t, double *ca, int nr, double k1, double k2, double k3, double kLoss, double *ct, double *cta, double *ctb)
Definition simulate.c:500
int simC3DIvs(double *t, double *ca1, double *ca2, double *cb, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double k1b, double k2b, double f, double vb, double fa, double *scpet, double *sct1, double *sct2, double *sct3, double *sct1b, double *sctab, double *sctvb)
Definition simulate.c:1387
int simC3s(double *t, double *ca, int nr, double k1, double k2, double k3, double k4, double k5, double k6, double *ct, double *cta, double *ctb, double *ctc)
Definition simulate.c:27