sim2cm.c

Go to the documentation of this file.

/*****************************************************************************/
#include "tpcclibConfig.h"
/*****************************************************************************/
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include <string.h>
/*****************************************************************************/
#include "tpccm.h"
/*****************************************************************************/
 
/*****************************************************************************/
int simC2(
  double *t,
  double *ca,
  const int nr,
  const double k1,
  const double k2,
  const double k3,
  const double k4,
  double *ct,
  double *cta,
  double *ctb
) {
  int i;
  double dt2, r, u, v;
  double cai, ca_last, t_last;
  double ct1, ct1_last, ct2, ct2_last;
  double ct1i, ct1i_last, ct2i, ct2i_last;
 
 
  /* Check for data */
  if(nr<2) return 1;
  if(t==NULL || ca==NULL || ct==NULL) return 2;
 
  /* Check parameters */
  if(k1<0.0) return 3;
 
  /* Calculate curves */
  t_last=0.0; if(t[0]<t_last) t_last=t[0];
  cai=ca_last=0.0;
  ct1_last=ct2_last=ct1i_last=ct2i_last=0.0;
  ct1=ct2=ct1i=ct2i=0.0;
  for(i=0; i<nr; i++) {
    /* delta time / 2 */
    dt2=0.5*(t[i]-t_last);
    /* calculate values */
    if(dt2<0.0) {
      return 5;
    } else if(dt2>0.0) {
      /* arterial integral */
      cai+=(ca[i]+ca_last)*dt2;
      /* Calculate partial results */
      r=1.0+k4*dt2;
      u=ct1i_last+dt2*ct1_last;
      v=ct2i_last+dt2*ct2_last;
      /* 1st tissue compartment and its integral */
      ct1 = ( k1*cai - (k2 + (k3/r))*u + (k4/r)*v ) / ( 1.0 + dt2*(k2 + (k3/r)) );
      ct1i = ct1i_last + dt2*(ct1_last+ct1);
      /* 2nd tissue compartment and its integral */
      ct2 = (k3*ct1i - k4*v) / r;
      ct2i = ct2i_last + dt2*(ct2_last+ct2);
    }
    /* copy values to argument arrays */
    ct[i]=ct1+ct2;
    if(cta!=NULL) cta[i]=ct1;
    if(ctb!=NULL) ctb[i]=ct2;
    /* prepare to the next loop */
    t_last=t[i]; ca_last=ca[i];
    ct1_last=ct1; ct1i_last=ct1i;
    ct2_last=ct2; ct2i_last=ct2i;
  }
 
  return 0;
}
/*****************************************************************************/
 
/*****************************************************************************/
int simC2_i(
  double *t,
  double *cai,
  const int nr,
  const double k1,
  const double k2,
  const double k3,
  const double k4,
  double *ct,
  double *cta,
  double *ctb
) {
  int i;
  double dt2, r, u, v;
  double t_last;
  double ct1, ct1_last, ct2, ct2_last;
  double ct1i, ct1i_last, ct2i, ct2i_last;
 
 
  /* Check for data */
  if(nr<2) return 1;
  if(t==NULL || cai==NULL || ct==NULL) return 2;
 
  /* Check parameters */
  if(k1<0.0) return 3;
 
  /* Calculate curves */
  t_last=0.0; if(t[0]<t_last) t_last=t[0];
  ct1_last=ct2_last=ct1i_last=ct2i_last=0.0;
  ct1=ct2=ct1i=ct2i=0.0;
  for(i=0; i<nr; i++) {
    /* delta time / 2 */
    dt2=0.5*(t[i]-t_last);
    /* calculate values */
    if(dt2<0.0) {
      return 5;
    } else if(dt2>0.0) {
      /* Calculate partial results */
      r=1.0+k4*dt2;
      u=ct1i_last+dt2*ct1_last;
      v=ct2i_last+dt2*ct2_last;
      /* 1st tissue compartment and its integral */
      ct1 = ( k1*cai[i] - (k2 + (k3/r))*u + (k4/r)*v ) / ( 1.0 + dt2*(k2 + (k3/r)) );
      ct1i = ct1i_last + dt2*(ct1_last+ct1);
      /* 2nd tissue compartment and its integral */
      ct2 = (k3*ct1i - k4*v) / r;
      ct2i = ct2i_last + dt2*(ct2_last+ct2);
    }
    /* copy values to argument arrays */
    ct[i]=ct1+ct2;
    if(cta!=NULL) cta[i]=ct1;
    if(ctb!=NULL) ctb[i]=ct2;
    /* prepare to the next loop */
    t_last=t[i];
    ct1_last=ct1; ct1i_last=ct1i;
    ct2_last=ct2; ct2i_last=ct2i;
  }
 
  return 0;
}
/*****************************************************************************/
 
/*****************************************************************************/