{.U+}

program model;

{This program models the operation of a Shiley M-2000 membrane lung during low
 blood flow extracorpereal CO2 removal.

 C. S. Tritt
 April 15, 1986}

{$I bldcalc2.pas}
{$I iolib.pas}

const

  phin: real = 7.42;
  temp: real = 37.0; {Degree C}
  bldnpco2: real = 42.0; {mmHg}
  bldnpo2: real = 160.0; {mmHg}
  hct: real = 40.0; {volume %}
  hbgm100: real = 15.0; {gm / 100 ml}
  baropres: real = 760; {mmHg}
  pno2: real = 21.0; {mole %}
  pnco2: real = 0.04; {mole %}
  gasflwin: real = 6.0; {l / min}
  bldflow: real = 0.5; {l / min}
  shuntco2 = 0.684; {Shunt based on CO2 for first flux approximation}
  ufixed = 5.26; {Fixed CO2 perm - ml / min mmHg}
  uslope = 2.63; {Slope of CO2 perm vs. BF line - ml / min mmHg l}
  uint = 2.99; {Intercept of CO2 perm vs. BF line - ml / min}
  gasresis = 1.15; {ML gas flow resistance - mmHg min / l}
  bldresis = 24.2; {ML blood flow resistance - mmHg min / l}
  pfors1 = 546.1; {Minimum PO2 for O2 shunt to equal 1 - mmHg}
  fo2slope = 0.001762; {Slope of O2 shunt vs. inlet PO2 line - / mmHg}
  fo2int = 0.038; {Intercept of O2 shunt vs. inlet PO2 line - unitless}

  fixed: boolean = false; {Fixed CO2 perm method used?}
  co2meth = 2; {1=Lutz's, 2=Kelman's or 3=Loeppky's dissociation curve used?}

var

  gasflwout, bldepo2, bldepco2, bldopo2, bldopco2, gasnpo2, gasnpco2, gasopo2,
    gasopco2, shunto2, flxo2, flxco2old, flxco2new, ph2o, qnco2, cnco2in,
    cno2in, cnco2eq, cno2eq, cnco2out, cno2out, permco2, be, pheq, phout,
    barcor, delpo2, delpco2, qno2, h2ofact, ln10, erro2, errco2, xno2, xnco2,
    pcorfact, gasnpres, gasdelp, blddelp, slope, intercep, flxinc, delold,
    delnew, delplm, delpbldo, delpbldn, gasopco2o, bldopco2o, d: real;

  goodio, masterr, error, phmeth2: boolean;

  loops: integer;

  species: strof14;

  entstrng: string[27];

  continue, ritefile, tryagain, phchar: char;

  parm: parmlist;

  desttag: text;

procedure display_heading;

begin
  clrscr;
  writeln('Chuck''s EXCCO2R Model for the Shiley M-2000');
  writeln;
  writeln('This program models the operation of the Shiley M-2000 membrane lung during');
  writeln('low blood flow extracorpereal CO2 removal.  It is based on experimental data');
  writeln('for pigs, but should also provide good results for other species.  Gas flows');
  writeln('should be between 0.5 and 12.0 l/min, blood flows between 0.25 and 1.00');
  writeln('l/min and gas outlet pressures (generally the barometric pressure) between');
  writeln('600 and 800 mmHg for best results.');
  writeln;
  writeln('To use this program answer the questions displayed on the screen.  To retain');
  writeln('a value just enter a RETURN in responce to the prompt.  Constants used');
  writeln('in the calculation of hemoglobin saturation should be in an ASCII text file.');
  writeln('Such data files have been provided for several common species.  If disk');
  writeln('output is selected, the output file will have the extension ''.out''.');
  writeln;
  writeln('It is the responsiblity of the user to determine if this program is suitable');
  writeln('for a particular application.  This program has been placed in the public');
  writeln('domain by the author and may be copied.  The financial support of the');
  writeln('Department of Anesthesology, The Ohio State University is acknowledged.');
  writeln;
  writeln('Version 3.0'); writeln('By C. S. Tritt'); writeln('April 14, 1986');
  writeln;
end;

procedure display_calculating;

{Ask user to wait 1 - 2 minutes for calculations}

var

  i: integer;

begin
  for i := 1 to 11 do writeln;
  write('Calculating.  Please wait.':51);
end;


begin {Mainline}

  entstrng := 'Enter new value or RETURN: ';
  ln10 := ln(10.0);
  delpo2 := 50.0; {mmHg}
  delpco2 := 10.0; {mmHg}

  display_heading;
  wait;
  clrscr;

  getlutz(species, parm);
  writeln;
  repeat
    write('Use linear pH-ln(Pco2) relationship? (y or n) ');
    read(kbd, phchar);
    phchar := upcase(phchar);
    if phchar = 'Y' then phmeth2 := false else phmeth2 := true;
  until (phchar = 'Y') or (phchar = 'N');
  writeln;
  getroute(ritefile, desttag);
  wait;
  clrscr;

  {Get conditions from user}
  repeat

    {Get global ML operating conditions}

    clrscr;
    writeln;
    putreal(baropres, 4, 0, 'gas outlet pressure in mmHg', con);
    getreal(baropres, entstrng);
    echoreal(baropres, 4, 0);
    writeln;
    putreal(hct, 4, 1, 'hematocrit in percent', con);
    getreal(hct, entstrng);
    echoreal(hct, 4, 1);
    writeln;
    putreal(hbgm100, 5, 2, 'hemoglobin in gm/100 ml', con);
    getreal(hbgm100, entstrng);
    echoreal(hbgm100, 5, 2);
    writeln;
    putreal(temp, 4, 1, 'temperature in C', con);
    getreal(temp, entstrng);
    echoreal(temp, 4, 1);
    writeln;
    clrscr;
    repeat

      {Get accessible ML operating conditions from user}

      writeln;
      writeln('Inlet gas conditions:');
      writeln;
      putreal(gasflwin, 5, 1, 'flow in l/min @ STPD', con);
      getreal(gasflwin, entstrng);
      echoreal(gasflwin, 5, 1);
      writeln;
      repeat
        putreal(pno2, 6, 1, 'oxygen mole percent', con);
        getreal(pno2, entstrng);
        echoreal(pno2, 6, 1);
        writeln;
        putreal(pnco2, 8, 3, 'carbon dioxide mole percent', con);
        getreal(pnco2, entstrng);
        echoreal(pnco2, 8, 3);
        writeln;
      until (pno2 + pnco2) <= 100.0;
      xno2 := pno2 / 100.0;
      xnco2 := pnco2 / 100.0;
      wait;
      clrscr;
      writeln;
      writeln('Intet blood conditions:');
      writeln;
      putreal(bldflow, 5, 2, 'blood flow in l/min', con);
      getreal(bldflow, entstrng);
      echoreal(bldflow, 5, 2);
      writeln;
      putreal(bldnpo2, 5, 1, 'PO2 in mmHg', con);
      getreal(bldnpo2, entstrng);
      echoreal(bldnpo2, 5, 1);
      writeln;
      putreal(bldnpco2, 6, 2, 'PCO2 in mmHg', con);
      getreal(bldnpco2, entstrng);
      echoreal(bldnpco2, 6, 2);
      writeln;
      putreal(phin, 5, 3, 'pH', con);
      getreal(phin, entstrng);
      echoreal(phin, 5, 3);
      writeln;
      wait;

      {Do preliminary calculations}

      clrscr;
      display_calculating;
      masterr := false;
      if fixed then permco2 := ufixed else permco2 := uslope * bldflow + uint;
      gasdelp := gasresis * gasflwin;
      blddelp := bldresis * bldflow;
      ph2o := fph2o(temp);
      barcor := baropres - ph2o;
      pcorfact := 760.0 * (273.0 + temp) / baropres / 273.0;
      gasnpres := baropres + gasdelp;
      h2ofact := 1.0 / (1.0 - ph2o / baropres);
      qno2 := xno2 * gasflwin;
      qnco2 := xnco2 * gasflwin;
      gasnpo2 := fgasxtop(pno2, gasnpres);
      gasnpco2 := fgasxtop(pnco2, gasnpres);
      gasopco2 := gasnpco2;

      if bldnpo2 > pfors1 then shunto2 := 1.0
      else shunto2 := fo2slope * bldnpo2 - fo2int;

      bldcon(cno2in, cnco2in, be, bldnpo2, bldnpco2, phin, temp, hct, hbgm100,
        hbgm100, co2meth, parm, error);
      masterr := masterr or error;
      bldepo2 := gasnpo2;
      bldopo2 := gasnpo2;
      bldepco2 := gasnpco2;
      phpco2si(slope, intercep, be, error);
      masterr := masterr or error;
      pheq := fphcalc2(bldepco2, slope, intercep, be, temp, phmeth2, error);
      masterr := masterr or error;
      phout := pheq + shuntco2 * (phin - pheq);
      bldcon(cno2eq, cnco2eq, d, bldepo2, bldepco2, pheq, temp, hct, hbgm100,
        hbgm100, co2meth, parm, error);
      masterr := masterr or error;
      cno2out := cno2eq + shunto2 * (cno2in - cno2eq);
      flxo2 := bldflow * (cno2out - cno2in) / 100.0;
      flxco2new := bldflow * (cnco2in - cnco2eq) / 100.0 * ( 1.0 - shuntco2);
      cnco2out := fconout(cnco2in, bldflow, flxco2new);
      bldopco2 := cnco2out / (2.22 * 0.03 * ( 1.0 + exp(ln10 * (phout -
        6.1))));

      flxinc := 0.050; {Set initial flux step in l/min}
      flxco2new :=  flxco2new - flxinc;
      delnew := - 1.0;

      loops := 0;

      {Start of interative routine}
      repeat
        loops := loops + 1;
        delold := delnew;
        gasopco2o := gasopco2;
        bldopco2o := bldopco2;

        if delnew <= 0.0 then flxco2new := flxco2new + flxinc
        else flxco2new := flxco2new - flxinc;
        if flxco2new < 0.0 then flxco2new := 0.0;

        cnco2out := fconout(cnco2in, bldflow, flxco2new);
        contops(bldopo2, bldopco2, phout, erro2, errco2, delpo2, delpco2, be,
          temp, hbgm100, hct, cno2out, cnco2out, parm, phmeth2, error);
        masterr := masterr or error;
        outgas(gasflwout, gasopco2, gasopo2, gasflwin, flxco2new, flxo2,
          pcorfact, h2ofact, baropres, qnco2, qno2, error);
        masterr := masterr or error;

        delpbldn := bldnpco2 - gasopco2;
        delpbldo := bldopco2 - gasnpco2;

        if delpbldn * delpbldo <= 0.0 then delnew := 1.0
        else
        begin
          delplm := (delpbldn - delpbldo) / ln(delpbldn / delpbldo);
          delnew := flxco2new - permco2 / 1000.0 * delplm;
        end;
        if delnew * delold < 0.0 then flxinc := flxinc / 2.0;

      until ((flxinc < 0.0005) and (abs(gasopco2 - gasopco2o)
        < 0.1) and (abs(bldopco2 - bldopco2o) < 0.1)) or (loops = 20);

      if loops = 20 then masterr := true;

      {Echo inputs}
      clrscr;
      writeln(desttag);
      putstring(species, 'Species: ', desttag);
      if phmeth2 then
        writeln(desttag, 'Thews'' pH-Pco2 relationship used.') else
        writeln(desttag, 'Linear pH-ln(Pco2) relationship used.');
      writeln(desttag);
      putreal(baropres, 4, 0, 'gas outlet pressure im mmHg', desttag);
      putreal(hct, 4, 1, 'haematocrit in percent', desttag);
      putreal(hbgm100, 5, 2, 'hemoglobin concentration in gm/100 ml', desttag);
      putreal(temp, 4, 1, 'temperature in C', desttag);
      writeln(desttag);
      writeln(desttag, 'Inlet conditions:');
      writeln(desttag);
      putreal(gasflwin, 5, 1, 'gas flow in l/min @ STPD', desttag);
      putreal(gasnpo2, 5, 1, 'gas PO2 in mmHg', desttag);
      putreal(gasnpco2, 6, 2, 'gas PCO2 in mmHg', desttag);
      writeln(desttag);
      putreal(bldflow, 5, 2, 'blood flow in l/min', desttag);
      putreal(bldnpo2, 5, 1, 'blood PO2 in mmHg', desttag);
      putreal(bldnpco2, 6, 2, 'blood PCO2 in mmHg', desttag);
      putreal(phin, 5, 3, 'pH', desttag);
      writeln(desttag);
      if ritefile = 'S' then wait;

      {Write results}
      clrscr;
      writeln(desttag, 'General results:');
      writeln(desttag);
      putreal(flxco2new, 6, 3, 'CO2 flux in l/min @ STPD (removed)', desttag);
      putreal(flxo2, 6, 3, 'O2 flux in l/min @ STPD (added)', desttag);
      writeln(desttag);
      putreal(gasdelp, 3, 0, 'gas pressure drop in mmHg', desttag);
      putreal(blddelp, 3, 0, 'blood pressure drop in mmHg', desttag);
      writeln(desttag);
      writeln(desttag, 'Outlet conditions:');
      writeln(desttag);
      putreal(gasflwout, 5, 1, 'gas flow in l/min @ BTPS', desttag);
      putreal(gasopo2, 5, 1, 'gas PO2 in mmHg', desttag);
      putreal(gasopco2, 6, 2, 'gas PCO2 in mmHg', desttag);
      writeln(desttag);
      putreal(bldopo2, 5, 1, 'blood PO2 in mmHg', desttag);
      putreal(bldopco2, 6, 2, 'blood PCO2 in mmHg', desttag);
      putreal(phout, 5, 3, 'pH', desttag);
      putreal(be, 7, 2, 'base excess in mEq/l', desttag);
      writeln(desttag);
      if ritefile = 'S' then wait;

      clrscr;
      {Get control information regarding repeat of calculations}
      writeln;
      write('Do you wish to continue? ( y or n ) ');
      read(kbd, continue);
      continue := upcase(continue);
      if continue <> 'N' then
      begin
        writeln;
        write('Do you wish change overall conditions? (y or n) ');
        read(kbd, tryagain); writeln;
        tryagain := upcase(tryagain);
      end
      else tryagain := 'N';
      clrscr;

    until (continue = 'N') or ((continue = 'Y') and (tryagain = 'Y'));
  until continue = 'N';

  {$I-}
  close(desttag);
  iocheck(goodio, 'Problem closing output file.');
  {$I+}
end.