TPCCLIB
Loading...
Searching...
No Matches
simkloss.c
Go to the documentation of this file.
1
4/*****************************************************************************/
5#include "tpcclibConfig.h"
6/*****************************************************************************/
7#include <stdio.h>
8#include <stdlib.h>
9#include <math.h>
10#include <time.h>
11#include <string.h>
12/*****************************************************************************/
13#include "tpccm.h"
14/*****************************************************************************/
15
16/*****************************************************************************/
36int simC3vpKLoss(
38 double *t,
40 double *ca,
42 double *cb,
44 const int nr,
46 const double k1,
48 const double k2,
50 const double k3,
52 const double k4,
54 const double k5,
56 const double k6,
58 const double kLoss,
60 const double f,
62 const double vb,
64 const double fa,
68 const int vvm,
70 double *cpet,
72 double *cta,
74 double *ctb,
76 double *ctc,
78 double *ctab,
80 double *ctvb
81) {
82 int i;
83 double dt2, b, c, d, w, z, u, va, vv;
84 double cai, ca_last, t_last, dct, cvb;
85 double ct1, ct1_last, ct2, ct2_last, ct3, ct3_last;
86 double ct1i, ct1i_last, ct2i, ct2i_last, ct3i, ct3i_last;
87
88
89 /* Check for data */
90 if(nr<2) return 1;
91 if(cpet==NULL) return 2;
92
93 /* Check parameters */
94 if(k1<0.0) return 3;
95 if(vb<0.0 || vb>=1.0) return 4;
96 if(fa<=0.0 || fa>1.0) return 5;
97 va=fa*vb; vv=(1.0-fa)*vb;
98
99 /* Calculate curves */
100 t_last=0.0; if(t[0]<t_last) t_last=t[0];
101 cai=ca_last=0.0;
102 ct1_last=ct2_last=ct3_last=ct1i_last=ct2i_last=ct3i_last=0.0;
103 ct1=ct2=ct3=ct1i=ct2i=ct3i=0.0;
104 for(i=0; i<nr; i++) {
105 /* delta time / 2 */
106 dt2=0.5*(t[i]-t_last);
107 /* calculate values */
108 if(dt2<0.0) {
109 return 5;
110 } else if(dt2>0.0) {
111 /* arterial integral */
112 cai+=(ca[i]+ca_last)*dt2;
113 /* Calculate partial results */
114 w = 1.0 + k4*dt2;
115 z = 1.0 + dt2*(k6 + kLoss);
116 u = k2 + k3 + k5 - k3*k4*dt2/w - k5*k6*dt2/z;
117 b = ct1i_last+dt2*ct1_last;
118 c = ct2i_last+dt2*ct2_last;
119 d = ct3i_last+dt2*ct3_last;
120 /* 1st tissue compartment and its integral */
121 ct1 = ( k1*cai - u*b + k4*c/w + k6*d/z ) / ( 1.0 + dt2*u);
122 ct1i = ct1i_last + dt2*(ct1_last+ct1);
123 /* 2nd tissue compartment and its integral */
124 ct2 = (k3*ct1i - k4*c) / w;
125 ct2i = ct2i_last + dt2*(ct2_last+ct2);
126 /* 3rd tissue compartment and its integral */
127 ct3 = (k5*ct1i - (k6 + kLoss)*d) / z;
128 ct3i = ct3i_last + dt2*(ct3_last+ct3);
129 }
130 /* Venous curve */
131 if(f>0.) {dct = k1*ca[i] - k2*ct1 - kLoss*ct3; cvb = cb[i] - dct/f;}
132 else cvb=cb[i];
133 /* copy values to argument arrays */
134 if(vvm==0) cpet[i]= va*cb[i] + vv*cvb + (1.0-vb)*(ct1+ct2+ct3);
135 else cpet[i]= va*cb[i] + vv*cvb + (ct1+ct2+ct3);
136 if(cta!=NULL) {cta[i]=ct1; if(vvm==0) cta[i]*=(1.0-vb);}
137 if(ctb!=NULL) {ctb[i]=ct2; if(vvm==0) ctb[i]*=(1.0-vb);}
138 if(ctc!=NULL) {ctc[i]=ct3; if(vvm==0) ctc[i]*=(1.0-vb);}
139 if(ctab!=NULL) {ctab[i]=va*cb[i];}
140 if(ctvb!=NULL) {ctvb[i]=vv*cvb;}
141 /* prepare to the next loop */
142 t_last=t[i]; ca_last=ca[i];
143 ct1_last=ct1; ct1i_last=ct1i;
144 ct2_last=ct2; ct2i_last=ct2i;
145 ct3_last=ct3; ct3i_last=ct3i;
146 }
147
148 return 0;
149}
150/*****************************************************************************/
151
152/*****************************************************************************/
169int simC2l(
171 double *t,
173 double *ca,
175 const int nr,
177 const double k1,
179 const double k2,
181 const double k3,
183 const double kLoss,
185 double *ct,
187 double *cta,
189 double *ctb
190) {
191 int i;
192 double b, c, dt2;
193 double cai, ca_last, t_last;
194 double ct1, ct1_last, ct2, ct2_last;
195 double ct1i, ct1i_last, ct2i, ct2i_last;
196
197
198 /* Check for data */
199 if(nr<2) return 1;
200 if(ct==NULL) return 2;
201
202 /* Check actual parameter number */
203 if(k1<0.0) return 3;
204
205 /* Calculate curves */
206 t_last=0.0; if(t[0]<t_last) t_last=t[0];
207 cai=ca_last=0.0;
208 ct1_last=ct2_last=ct1i_last=ct2i_last=ct1=ct2=ct1i=ct2i=0.0;
209 for(i=0; i<nr; i++) {
210 /* delta time / 2 */
211 dt2=0.5*(t[i]-t_last);
212 /* calculate values */
213 if(dt2<0.0) {
214 return 5;
215 } else if(dt2>0.0) {
216 /* arterial integral */
217 cai+=(ca[i]+ca_last)*dt2;
218 /* partial results */
219 b=ct1i_last+dt2*ct1_last;
220 c=ct2i_last+dt2*ct2_last;
221 /* 1st tissue compartment and its integral */
222 ct1 = (k1*cai - (k2+k3)*b) / (1.0 + (k2+k3)*dt2 );
223 ct1i = ct1i_last + dt2*(ct1_last+ct1);
224 /* 2nd tissue compartment and its integral */
225 ct2 = (k3*ct1i - kLoss*c) / (1.0 + kLoss*dt2);
226 ct2i = ct2i_last + dt2*(ct2_last+ct2);
227 }
228 /* copy values to argument arrays */
229 ct[i]=ct1+ct2;
230 if(cta!=NULL) cta[i]=ct1;
231 if(ctb!=NULL) ctb[i]=ct2;
232 /* prepare to the next loop */
233 t_last=t[i]; ca_last=ca[i];
234 ct1_last=ct1; ct1i_last=ct1i;
235 ct2_last=ct2; ct2i_last=ct2i;
236 }
237
238 return 0;
239}
240/*****************************************************************************/
241
242/*****************************************************************************/
261int simC2vl(
263 double *t,
265 double *ca,
267 double *cb,
269 const int nr,
271 const double k1,
273 const double k2,
275 const double k3,
277 const double kL,
279 const double f,
281 const double vb,
283 const double fa,
287 const int vvm,
289 double *cpet,
291 double *cta,
293 double *ctb,
295 double *ctab,
297 double *ctvb
298) {
299 int i;
300 double dt2, b, c, va, vv;
301 double cai, ca_last, t_last, dct, cvb;
302 double ct1, ct1_last, ct2, ct2_last;
303 double ct1i, ct1i_last, ct2i, ct2i_last;
304
305
306 /* Check for data */
307 if(nr<2) return 1;
308 if(cpet==NULL) return 2;
309
310 /* Check parameters */
311 if(k1<0.0) return 3;
312 if(vb<0.0 || vb>=1.0) return 4;
313 if(fa<=0.0 || fa>1.0) return 5;
314 va=fa*vb; vv=(1.0-fa)*vb;
315
316 /* Calculate curves */
317 t_last=0.0; if(t[0]<t_last) t_last=t[0];
318 cai=ca_last=0.0;
319 ct1_last=ct2_last=ct1i_last=ct2i_last=ct1=ct2=ct1i=ct2i=0.0;
320 for(i=0; i<nr; i++) {
321 /* delta time / 2 */
322 dt2=0.5*(t[i]-t_last);
323 /* calculate values */
324 if(dt2<0.0) {
325 return 5;
326 } else if(dt2>0.0) {
327 /* arterial integral */
328 cai+=(ca[i]+ca_last)*dt2;
329 /* Calculate partial results */
330 b=ct1i_last+dt2*ct1_last;
331 c=ct2i_last+dt2*ct2_last;
332 /* 1st tissue compartment and its integral */
333 ct1 = (k1*cai - (k2+k3)*b) / (1.0 + (k2+k3)*dt2 );
334 ct1i = ct1i_last + dt2*(ct1_last+ct1);
335 /* 2nd tissue compartment and its integral */
336 ct2 = (k3*ct1i - kL*c) / (1.0 + kL*dt2);
337 ct2i = ct2i_last + dt2*(ct2_last+ct2);
338 }
339 /* Venous curve */
340 if(f>0.) {dct = k1*ca[i] - k2*ct1 - kL*ct2; cvb = cb[i] - dct/f;}
341 else cvb=cb[i];
342
343 /* copy values to argument arrays */
344 if(vvm==0) cpet[i]= va*cb[i] + vv*cvb + (1.0-vb)*(ct1+ct2);
345 else cpet[i]= va*cb[i] + vv*cvb + (ct1+ct2);
346 if(cta!=NULL) {cta[i]=ct1; if(vvm==0) cta[i]*=(1.0-vb);}
347 if(ctb!=NULL) {ctb[i]=ct2; if(vvm==0) ctb[i]*=(1.0-vb);}
348 if(ctab!=NULL) {ctab[i]=va*cb[i];}
349 if(ctvb!=NULL) {ctvb[i]=vv*cvb;}
350 /* prepare to the next loop */
351 t_last=t[i]; ca_last=ca[i];
352 ct1_last=ct1; ct1i_last=ct1i;
353 ct2_last=ct2; ct2i_last=ct2i;
354 }
355
356 return 0;
357}
358/*****************************************************************************/
359
360/*****************************************************************************/
simC3vpKLoss
int simC3vpKLoss(double *t, double *ca, double *cb, const int nr, const double k1, const double k2, const double k3, const double k4, const double k5, const double k6, const double kLoss, const double f, const double vb, const double fa, const int vvm, double *cpet, double *cta, double *ctb, double *ctc, double *ctab, double *ctvb)
Definition simkloss.c:36
simC2l
int simC2l(double *t, double *ca, const int nr, const double k1, const double k2, const double k3, const double kLoss, double *ct, double *cta, double *ctb)
Definition simkloss.c:169
simC2vl
int simC2vl(double *t, double *ca, double *cb, const int nr, const double k1, const double k2, const double k3, const double kL, const double f, const double vb, const double fa, const int vvm, double *cpet, double *cta, double *ctb, double *ctab, double *ctvb)
Definition simkloss.c:261
tpccm.h
Header file for libtpccm.
Generated on for TPCCLIB by doxygen 1.15.0