unit Unit2;

interface

uses
  Winapi.Windows, Winapi.Messages, System.SysUtils, System.Variants, System.Classes, Vcl.Graphics,
  Vcl.Controls, Vcl.Forms, Vcl.Dialogs, Vcl.StdCtrls, Vcl.Samples.Spin,
  Vcl.ExtCtrls, Vcl.Menus, IOUtils;

const
unknown = -9999;
tab = #9;
maxNumber = 1000000000000000;
maxAge = 100;
minWedAge = 15.5;
maxduration = 500;
FertilityThreshold = 45;
MaxNumberOfData = 24;
NumberOfPopData = 24;
NumberOfHousesData = 6;
PopFit = 1;
HousesFit = 2;

PopEngland: array[1..NumberOfPopData,1..2] of real =
((1541, 2.830), (1546, 2.908), (1551, 3.065), (1556, 3.213),
(1561, 3.036), (1566, 3.174), (1571, 3.310), (1576, 3.448),
(1581, 3.631), (1586, 3.841), (1591, 3.938), (1596, 4.057),
(1601, 4.162), (1606, 4.310), (1611, 4.476), (1616, 4.568),
(1621, 4.745), (1626, 4.762), (1631, 4.926), (1636, 5.090),
(1641, 5.130), (1646, 5.231), (1651, 5.308), (1656, 5.391));

HousesEngland: array[1..NumberOfHousesData,1..3] of real =
((1116, 0.535, 76), (1181, 0.671, 67), (1226, 0.985, 49),
(1250, 0.985, 111), (1267, 1.152, 51), (1295, 1.077, 398));



type Distribution = array[0..maxAge] of real;

type
  TForm2 = class(TForm)
    Panel1: TPanel;
    Label1: TLabel;
    Label2: TLabel;
    Label3: TLabel;
    Label4: TLabel;
    Edit1: TEdit;
    Edit2: TEdit;
    Edit3: TEdit;
    Label7: TLabel;
    Panel2: TPanel;
    Label8: TLabel;
    Label9: TLabel;
    Label10: TLabel;
    Edit7: TEdit;
    Label11: TLabel;
    Button1: TButton;
    Label12: TLabel;
    Label13: TLabel;
    Label14: TLabel;
    Edit6: TEdit;
    Panel3: TPanel;
    Label5: TLabel;
    Label6: TLabel;
    Edit4: TEdit;
    Label15: TLabel;
    Edit5: TEdit;
    Label19: TLabel;
    SpinEdit4: TSpinEdit;
    Panel4: TPanel;
    Label16: TLabel;
    Memo1: TMemo;
    Memo2: TMemo;
    Label17: TLabel;
    Label18: TLabel;
    CheckBox1: TCheckBox;
    SpinEdit1: TSpinEdit;
    SpinEdit2: TSpinEdit;
    Label21: TLabel;
    Panel5: TPanel;
    Label20: TLabel;
    Label22: TLabel;
    Label23: TLabel;
    Edit8: TEdit;
    Edit9: TEdit;
    Button2: TButton;
    Button3: TButton;
    Button4: TButton;
    SpinEdit3: TSpinEdit;
    Spinedit5: TSpinEdit;
    Label24: TLabel;
    Label25: TLabel;
    procedure FormCreate(Sender: TObject);
    procedure Button1Click(Sender: TObject);
    procedure Button2Click(Sender: TObject);
    procedure Button3Click(Sender: TObject);
    procedure Button4Click(Sender: TObject);
  private
    { Private declarations }
  public
    { Public declarations }
    logdir:string;
    logfile:string;
    //parameters
    TypeOfFit:integer;
    beta:real;
    theta:real;
    ohmega:real;
    BlackDeathMortality:real;
    gamma:real;
    WeddingMean:integer;
    WeddingK :real;
    BpY:real;
    kappa:real;
    labda:real;
    growthrateHouseholds:real;
    Growthtype: integer;
    Householdsize:real;
    Pbegin:real;
    //variables
    tbegin:integer;
    t:integer;
    duration:integer;
    Positions:real;
    PtotalJan:real;
    PtotalDec:real;
    Pt:distribution;
    HouseholdsTotal:real;
    Households: distribution;
    BachelorsTotal:real;
    Bachelors:distribution;
     Weddings:distribution;
    nQx :distribution;
    U14t:real;
    gt:real;
    births :real;
    deaths:real;
    deathsHouseholds:real;
    MeanAgeMarried:real;
    MeanAgeWedded:real;
    NumberOfData:integer;
    error:real;
    Ngrowthrate:integer;
    Npopulation:integer;
    function WriteReal(r:real;decimalen:integer):string;
    procedure InitiateVariables;
    procedure InitiateParameters;
    function Survivorship(x: Real):real;
    procedure BuildPreCrisisMortality;
    procedure BuildPopulationJanuaryStartyear;
    procedure BuildPopulationJanuary;
    Procedure BuildStartHouseholds;
    procedure CalculateTotal(var Dtotal:real;D:distribution);
    procedure CalculateU14;
    procedure CalculateBirths;
    procedure CalculateDeaths(var D:distribution; var rdeaths, rtotal:real);
    procedure CalculateBachelors;
    procedure CalculateWeddings;
    procedure CalculatePositions;
    procedure CalculateLossPositions;
    procedure CalculateGrowthrateHouseholds;
    procedure CalculateMeanAgeMarried;
    procedure CalculateMeanAgeWedded;
    procedure SimulationAdaptingBirths;
    procedure TestSimulationAdaptingBirths;
    function CDF(x:real):real;
    function fit:real;
    {function GaussCurve(i:integer):real; }
    function CoaleMcNeil(x:real):real;
    procedure Data(Px,tx:real);
    procedure HousesData(Px,tx:real);
    procedure Proxy(Px,tx:real);
  end;

var
  Form2: TForm2;
  TestYears: array[1..MaxNumberOfData] of integer;
  TestArray: array[1..MaxNumberOfData,0..2] of real;


implementation

{$R *.dfm}

uses Unit1, NelderMead;


procedure TForm2.Button2Click(Sender: TObject);
var i, iteration: integer;
    S: TSimplex;
    ss:string;
begin
  typeOfFit := PopFit;
  NumberOfData := NumberOfPopData;
  memo1.Lines.Clear;
  InitiateVariables;
  InitiateParameters;
  form1.InitiateGraphs;
  for i:= 1 to NumberOfData
  do begin
     TestYears[i] := Round(PopEngland[i,1]);
     TestArray[i,0] := 1; // weights
     TestArray[i,1] := PopEngland[i,2];
     end;

  S := TSimplex.New(kappa*1000,labda*100,PBegin,Householdsize,3);
  S.StartSimplex;
  Iteration := 0;
  While ((iteration < 5000) and (S.Margin > dSimplex/10000))
  do begin
     S.Iteration;
     Inc(iteration,1);
     ss := 'iteration= ' + writereal(iteration,0) + ' error= ' + writereal(S.CentreValue,2) +  ' k= ' + writereal(kappa,5) + ' l= ' + writereal(labda,4) + ' Pbegin= ' + writereal(PBegin,2);
     caption := ss;
     memo1.Lines.Add(ss);
     end;
  S.Free; S:= nil;

  form1.TemporarySave;
  form1.Caption := caption;
  form1.Visible := true;
end;

procedure TForm2.Button3Click(Sender: TObject);
var i1,i2,i3,i: integer;
    kappaStart, labdaStart, PBeginStart, Bestfit, r:real;
    ss:string;
begin
  TypeOfFit := PopFit;
  NumberOfData := NumberOfPopdata;
  memo1.Lines.Clear;
  InitiateVariables;
  InitiateParameters;
  form1.InitiateGraphs;
  for i:= 1 to NumberOfData
  do begin
     TestYears[i] := Round(PopEngland[i,1]);
     TestArray[i,0] := 1; // Weights
     TestArray[i,1] := PopEngland[i,2];
     end;
 kappaStart := kappa;
 labdaStart := labda;
 PBeginStart := PBegin;
 BestFit := maxNumber;
 For i1 := -Ngrowthrate to Ngrowthrate
 do begin
    {BestFit := maxNumber;}
    for i2 := -Ngrowthrate to Ngrowthrate
    do for i3 := -Npopulation to Npopulation
       do begin
          kappa := kappaStart + i1 * 0.00001;
          labda := labdaStart + i2 * 0.0001;
          PBegin := PBeginStart + i3 * 0.01;
          PtotalJan := PBegin;
          CalculateGrowthrateHouseholds;
          Householdsize := 5.87 - 0.821 * Growthratehouseholds*100;
          TestSimulationAdaptingBirths;
          r := fit;
          ss :=  ' error= ' + writereal(r,5) + tab +  'k= ' + writereal(kappa,5) + tab + ' l=' + writereal(labda,4) + tab + 'Pbegin= ' + writereal(PBegin,2);
          caption := ss;
          If r < Bestfit
          then begin
               memo1.Lines.Add(ss);
               Bestfit := r;
               end;
          end;
    end;
  form1.TemporarySave;
  form1.Caption := caption;
  form1.Visible := true;
end;

procedure TForm2.Button4Click(Sender: TObject);
var i1,i2,i3,i: integer;
    kappaStart, labdaStart, PBeginStart, Bestfit, r:real;
    ss:string;
begin
  TypeOfFit := HousesFit;
  NumberOfData := NumberOfHousesData;
  memo1.Lines.Clear;
  InitiateVariables;
  InitiateParameters;
  form1.InitiateGraphs;
  for i:= 1 to NumberOfData
  do begin
     TestYears[i] := Round(HousesEngland[i,1]);
     TestArray[i,0] := HousesEngland[i,3]; //weights
     TestArray[i,1] := HousesEngland[i,2];
     end;
 kappaStart := kappa;
 labdaStart := labda;
 PBeginStart := PBegin;
 BestFit := maxNumber;
 For i1 := -Ngrowthrate to Ngrowthrate
 do begin
    {BestFit := maxNumber;}
    for i2 := -Ngrowthrate to Ngrowthrate
    do for i3 := -Npopulation to Npopulation
       do begin
          kappa := kappaStart + i1 * 0.00001;
          labda := labdaStart + i2 * 0.0001;
          PBegin := PBeginStart + i3 * 0.01;
          PtotalJan := PBegin;
          CalculateGrowthrateHouseholds;
          Householdsize := 5.87 - 0.821 * Growthratehouseholds*100;
          TestSimulationAdaptingBirths;
          r := fit;
          ss :=  ' error= ' + writereal(r,5) + tab +  'k= ' + writereal(kappa,5) + tab +  'l= ' + writereal(labda,4) + tab + 'Pbegin= ' + writereal(PBegin,2);
          caption := ss;
          If r < Bestfit
          then begin
               memo1.Lines.Add(ss);
               Bestfit := r;
               end;
          end;
    end;
  form1.TemporarySave;
  form1.Caption := caption;
  form1.Visible := true;

end;

function Tform2.fit:real;
var i:integer;
    r, rtotal, weightsTotal :real;
begin
  rtotal:= 0;
  weightsTotal := 0;
  for i:= 1 to NumberOfData
  do begin
     r := sqr(TestArray[i,1] - TestArray[i,2]);
     rtotal := rtotal + TestArray[i,0] * r;
     weightsTotal := weightsTotal  + TestArray[i,0];
     end;
  result := SQRT(rtotal/weightsTotal);
end;

procedure Tform2.TestSimulationAdaptingBirths;
var PoldJan , Px, tx:real;
    i,j:integer;
    s:string;
begin
  t := tbegin;
  BuildPopulationJanuaryStartyear;
  CalculateDeaths(Pt,deaths,Ptotaldec);
  CalculatePositions;
  CalculateU14;
  {gt := growthrateHouseholds;}
  try
    for i:= 1 to duration
    do begin
       if ((ABS(PtotalJan) < maxNumber) and (ABS(U14t) < maxNumber))
       then begin
            t := t + 1;
            PoldJan := PtotalJan;
            BuildPopulationJanuary;
            CalculateDeaths(Pt,deaths,Ptotaldec);
            CalculateWeddings;
            CalculatePositions;
            CalculateTotal(HouseholdsTotal,Households);
            gt := (PtotalJan/PoldJan -1) * 100;
            CalculateU14;
            CalculateMeanAgeMarried;
            CalculateMeanAgeWedded;
            for j:= 1 to NumberOfData
            do if Testyears[j] = t
            then begin
                 case TypeOfFit of
                  PopFit :  TestArray[j,2] := Ptotaljan;
                  HousesFit: TestArray[j,2] := Ptotaljan / Householdsize;
                  end;
                 end;
            end;
       end;
  except
   form2.caption := inttostr(i);
  end;

end;



{function TForm2.GaussCurve(i:integer):real;
begin
   result := exp(-0.5 *sqr(i-WeddingTop)/sqr(WeddingSD));
end;}

procedure Tform2.CalculateGrowthrateHouseholds;
var wagereal:real;
begin
   GrowthrateHouseholds := -kappa * PtotalJan + labda;
end;

procedure TForm2.CalculateLossPositions;
begin
end;


procedure Tform2.CalculatePositions;
var vacancies:real;
begin
  vacancies := Positions - HouseholdsTotal;
  CalculateGrowthrateHouseholds;
  If vacancies < 0.001 * Positions 
  then Positions := Positions * (1 + growthrateHouseholds);
end;


procedure Tform2.CalculateTotal(var Dtotal:real; D:distribution);
var i:integer;
begin
  Dtotal := 0;
  for i:= 0 to maxAge do Dtotal := Dtotal + D[i];
end;



function Tform2.CoaleMcNeil(x:real):real;
var r:real;
begin
 r := {0.19465/WeddingK * par[1]*} exp(-0.174/WeddingK*(x-minWedAge-6.06*WeddingK)
                                         - exp(-0.288/WeddingK*(x-minWedAge-6.06*WeddingK)));
 result := r;
end;


function TForm2.CDF(x: real):real;
var r, rtot, phi:real;
    i:integer;
begin
  r := x;
  rtot := r;
  for i:= 1 to 50
  do begin
     r := r * x * x /(2*i+1);
     rtot := rtot + r;
     end;
  phi := 0.5 + 1/sqrt(2*pi) * exp(-sqr(x)/2) * rtot;
  result := phi;
end;


function TForm2.WriteReal(r:real;decimalen:integer):string;
var s:string;
begin
 str(r:2:decimalen,s);
 result := s;
end;




procedure Tform2.InitiateVariables;
var code:integer;
begin
  tBegin := spinedit1.value;
  duration := spinedit2.value;
  val(edit7.text,PBegin,code);
end;


procedure Tform2.InitiateParameters;
var code:integer;
begin
  val(edit1.Text,beta,code);
  val(edit2.Text,theta,code);
  ohmega := 119.3*beta  -1.01*theta  -14.5;
  If ohmega < 2.7 then ohmega := 2.7;
  val(edit3.Text,BlackDeathMortality,code);
  val(edit6.Text,gamma,code);
  val(edit4.Text,Householdsize,code);
  val(edit5.Text,BpY,code);
  val(edit8.Text,kappa,code);
  val(edit9.Text,labda,code);
  Ngrowthrate := spinedit3.Value;
  NPopulation := spinedit5.value;
  {WeddingTop := spinedit4.Value;
  WeddingSD := spinedit5.Value; }
  WeddingMean := spinedit4.Value;
  WeddingK := (WeddingMean - minWedAge)/11.36;
end;


function TForm2.Survivorship(x: Real):real;
var r0, r , rhulp:real;
    s:string;
begin
  //power formula for immature term, logistic for elderly term:
  r0 := 1 / exp( -1/(1-ohmega) * exp((1-ohmega)*ln(1.5))
                 - exp(theta)*ln(1 + exp(-16))
                 -1/beta * ln(1 + exp(-beta*100)));

  r := r0 * exp( -1/(1-ohmega) * exp((1-ohmega)*ln(x+1.5))
                 - exp(theta)*x
                 - exp(theta)*ln(1 + exp(x - 16))
                 -1/beta * ln(1 + exp(beta * (x - 100))));
  result := r;
end;

procedure TForm2.BuildPreCrisisMortality;
var i:integer;
    r1,r2:real;
begin
  for i:= 0 to maxAge -1
  do begin
     r1 := Survivorship(i);
     r2 := Survivorship(i+1);
     nQx[i] := (r1 - r2)/r1;
     end;
  nQx[maxAge] := 1;
end;


procedure TForm2.BuildPopulationJanuaryStartyear;
var i:integer;
    r:real;
begin
  PtotalJan := Pbegin;
  BuildPreCrisisMortality;
  CalculateGrowthrateHouseholds;
  for i:= 0 to maxAge do Pt[i] := survivorship(i) * exp(-growthrateHouseholds*i); //assume steady growth
  CalculateTotal(PtotalJan, Pt);
  for i:= 0 to maxAge do Pt[i] := Pt[i]/PtotalJan * Pbegin;
  births := Pt[0];
  CalculateTotal(PtotalJan, Pt);
  {for i:= 0 to maxAge do Pt[i] := Pt[i] * (1 - nQx[i]);}
  BuildstartHouseholds;
end;


procedure Tform2.BuildStartHouseholds;
var i:integer;
    Rtot, Ri : real;
    RR:array[0..maxAge] of real;
begin
  //calculate start positions
  Positions := Pbegin / Householdsize;

  HouseholdsTotal := Positions;
  for i:= 0 to maxAge do Households[i] := 0;

  for i:= 0 to round(minWedAge) -1 do RR[i] := 0;
  for i:= round(minWedAge) to MaxAge
  do RR[i] := (1-nQx[i])* RR[i-1] * exp(-growthrateHouseholds)
            + CoaleMcNeil(i)* exp(-growthrateHouseholds*(i-round(minWedAge))) {GaussCurve(i)};
  Rtot := 0;
  for i:= round(minWedAge) to MaxAge do Rtot := Rtot + RR[i];
  for i:= round(minWedAge) to MaxAge do Households[i] := RR[i]/ Rtot * HouseholdsTotal;
end;




procedure TForm2.BuildPopulationJanuary;
var i:integer;
begin
  for i:= maxAge downto 1 do Pt[i] := Pt[i-1];
  CalculateBirths;
  CalculateTotal(PtotalJan, Pt);
   Householdsize := PtotalJan/Householdstotal;
  for i:= maxAge downto 1 do Households[i] := Households[i-1];
end;



procedure TForm2.Button1Click(Sender: TObject);
begin
  InitiateVariables;
  InitiateParameters;
  form1.InitiateGraphs;
  SimulationAdaptingBirths;
  form1.TemporarySave;
  form1.Caption := caption;
  form1.Visible := true;
end;



procedure TForm2.CalculateU14;
var r:real;
    i:integer;
begin
  r:= 0;
  for i:= 0 to 13 do r := r + Pt[i];
  U14t := r/PtotalDec * 100;
end;


procedure TForm2.FormCreate(Sender: TObject);
var i:integer;
    s:string;
    F:textfile;
begin
  s := DateToStr(Date) + ' ' + Timetostr(Time);
  for i:= 1 to length(s) do if s[i] in ['/',':'] then s[i] := '_';
  logdir :=  getcurrentdir + '\logfiles\';
  logfile := logdir + 'log' + s +  '.txt';
  If not Tdirectory.exists(logdir) then Tdirectory.CreateDirectory(logdir);
  assignfile(F,logfile);
  rewrite(F);
  closefile(F);
end;


procedure TForm2.CalculateBirths;
var i,n:integer;
    rr:real;
begin
  births := 0;

  for i:= 15 to FertilityThreshold
  do begin
     rr := (1 - exp( 4.6 * ( (i-15) /(FertilityThreshold-15) -1)));
     births := births + rr * Households[i] * BpY;
     end;

  {for i:= 0 to FertDeclineAge -1 do births := births + Households[i] * (1-infertilePerc/100)/birthinterval;
  for i:= FertDeclineAge to FertilityThreshold
  do begin
     r := (fertilityThreshold - i) / (fertilityThreshold - FertDeclineAge);
     births := births + r * Households[i] * (1-infertilePerc/100)/birthinterval;
     end;}

  Pt[0] := births;
end;


 
procedure Tform2.CalculateWeddings;
var i:integer;
    vacancies, WeddingsTotal, rhulp, r, rr, rtot, C:real;
begin
  CalculateDeaths(Households,deathsHouseholds,Householdstotal);
  CalculateLossPositions;
  CalculateBachelors;
  vacancies := Positions - HouseholdsTotal;
  If ((vacancies > 0) and (BachelorsTotal > 0))
  then begin
       WeddingsTotal := vacancies;
       for i:= 0 to maxAge do Weddings[i] := 0;


       C:= 0;
       for i:= round(minWedAge) to maxAge
       do  begin
           r := CoaleMcNeil(i) * exp(-growthrateHouseholds*(i-round(minWedAge)));
           C := C + r;
           end;


       rtot := 0;
       for i:= round(minWedAge) to maxAge
       do begin
          r := CoaleMcNeil(i)/C *exp(-growthrateHouseholds*(i-round(minWedAge))) ;
          rr := r *  WeddingsTotal;
          If rr > Bachelors[i] then rr := Bachelors[i];
          rtot := rtot + rr ;
          Weddings[i] := rr;
          end;

       for i:= 0 to round(minWedAge)-1 do Households[i] := 0;
       for i:= round(minWedAge) to maxAge do Households[i] := Households[i] + Weddings[i];
       end;

  CalculateTotal(Householdstotal,Households);
  Calculatebachelors;
end;


procedure Tform2.CalculateBachelors;
var i:integer;
    r:real;
begin
  for i:= round(minWedAge) to maxAge
  do begin
     Bachelors[i] := Pt[i]/2 - Households[i];
     if Bachelors[i] < 0 then Bachelors[i] := 0;
     end;
  BachelorsTotal := 0;
  for i:= round(minWedAge) to maxAge do BachelorsTotal := BachelorsTotal + Bachelors[i];
end;



procedure TForm2.CalculateDeaths(var D:distribution; var rdeaths, rtotal:real);
var i:integer;
    r, mp:real;
begin
  rdeaths := 0;
  for i:= 0 to maxAge
  do begin
     case t of
       1317: r := 0.10 + nQx[i];
       1349, 1360, 1369, 1375, 1390, 1400, 1415, 1430, 1445
                              : begin
                                mp := BlackDeathMortality/100 * exp(-gamma * (t-1349));
                                if mp < 0.15 then mp := 0.15;
                                if i=0
                                then r:= 2 * mp + nQx[i]
                                else r:= mp + nQx[i];
                                end;
       else r:= nQx[i];
     end;
     if r > 1 then r:= 1;
     D[i] := D[i] * (1-r);
     rdeaths := rdeaths + r * D[i];
     end;
  CalculateTotal(rtotal,D);
end;


procedure Tform2.CalculateMeanAgeMarried;
var i:integer;
    r:real;
begin
  r := 0;
  for i:= 0 to MaxAge do r := r + i * Households[i];
  MeanAgeMarried := r / Householdstotal;
end;

procedure Tform2.CalculateMeanAgeWedded;
var i:integer;
    r, rr:real;
begin
  r := 0;
  rr := 0;
  for i:= 0 to MaxAge do r := r + i * Weddings[i];
  for i:= 0 to MaxAge do rr := rr + Weddings[i];
  MeanAgeWedded := r / rr;
end;


procedure Tform2.SimulationAdaptingBirths;
var PoldJan , Px, tx:real;
    i,j, n1,m1,n2,m2,n3,m3,n4,m4,n5,m5,n6,m6,n7,m7,n8,m8:integer;
    s:string;
begin
  n1 := unknown;  m1 := unknown;
  n2 := unknown;  m2 := unknown;
  n3 := unknown;  m3 := unknown;
  n4 := unknown;  m4 := unknown;
  n5 := unknown;  m5 := unknown;
  n6 := unknown;  m6 := unknown;
  n7 := unknown;  m7 := unknown;
  n8 := unknown;  m8 := unknown;
  memo1.Clear;
  memo1.lines.Add('year' + tab + 'Pop.' + tab + 'U 14' + tab + 'Births' + tab + 'Deaths' + tab + 'Deaths Hoseholds' +tab + 'growth(%)'
      + tab + 'Households' + tab + 'Mean Age Married' + tab + 'Mean Age Wedded' + tab + 'Household size'+ tab + 'Growthttype');
  memo2.Clear;
  memo2.Lines.Add('age' + Tab + '20' + tab + '25' + tab + '30' + tab + '35'
                         + tab + '40' + tab + '45' + tab + '50' + tab + '55'
                         + tab + '60' + tab + '65' + tab + '70' + tab + '75');
  t := tbegin;
  BuildPopulationJanuaryStartyear;
  CalculateDeaths(Pt,deaths,Ptotaldec);
  CalculatePositions;
  CalculateU14;
  try
    for i:= 1 to duration
    do begin
       if ((ABS(PtotalJan) < maxNumber) and (ABS(U14t) < maxNumber))
       then begin
            t := t + 1;
            PoldJan := PtotalJan;
            BuildPopulationJanuary;
            CalculateDeaths(Pt,deaths,Ptotaldec);
            CalculateWeddings;
            CalculatePositions;
            CalculateTotal(HouseholdsTotal,Households);
            gt := (PtotalJan/PoldJan -1) * 100;
            CalculateU14;
            CalculateMeanAgeMarried;
            CalculateMeanAgeWedded;
            memo1.lines.Add(inttostr(t) + tab + writereal(PtotalJan,3) + tab + writereal(U14t,0)
                  + tab + writereal(births,2) + tab + writereal(deaths,2) + tab + writereal(deathsHouseholds,3)
                  + tab +tab + writereal(gt,2) + tab + tab + writereal(HouseholdsTotal,3)
                  + tab + tab + writereal(MeanAgeMarried,0) + tab + tab + writereal(MeanAgeWedded,0)
                  + tab + tab + writereal(Householdsize,1) + tab + inttostr(Growthtype));
            s := '';
            for j:= 4 to 15
            do if Pt[j*5] > 0
               then s := s + tab + writereal(100* Bachelors[j*5]/(Pt[j*5]/2),0)
               else s := s + tab + '-';
            memo2.Lines.Add(inttostr(t) + s);
            If t > tBegin
            then begin
                 form1.DrawpointPt(PtotalJan,t,n2,m2,clblack,1,true);
                 form1.DrawpointU14(U14t,t,n3,m3,clblack,1,true);
                 form1.DrawpointFt(births,t,n4,m4,clblack,1,true);
                 form1.DrawpointGt(gt,t,n5,m5,clblack,1,true);
                 form1.DrawpointHt(HouseholdsTotal,t,n6,m6,clblack,1,true);
                 form1.DrawPointAt(MeanAgeMarried,t,n7,m7,clblack,1,true);
                 form1.DrawpointHSt(Householdsize,t,n8,m8,clblack,1,true);
                 end;
            end;
       end;
  except
   form2.caption := inttostr(i);
  end;
end;



end.




unit Unit1;

interface

uses
  Winapi.Windows, Winapi.Messages, System.SysUtils, System.Variants, System.Classes, Vcl.Graphics,
  Vcl.Controls, Vcl.Forms, Vcl.Dialogs, Vcl.ExtCtrls, Vcl.ComCtrls, Vcl.Menus;

const
small = 0.0000001;
hheight = 400;
wwidth = 900;

type
  TForm1 = class(TForm)
    MainMenu1: TMainMenu;
    Panel1: TPanel;
    PageControl1: TPageControl;
    TabSheet1: TTabSheet;
    ImagePt: TImage;
    TabSheet2: TTabSheet;
    ImageU14: TImage;
    TabSheet3: TTabSheet;
    ImageFt: TImage;
    TabSheet4: TTabSheet;
    Imagegt: TImage;
    TabSheet5: TTabSheet;
    ImageHt: TImage;
    file1: TMenuItem;
    savegraph1: TMenuItem;
    TabSheet6: TTabSheet;
    TabSheet7: TTabSheet;
    ImageAt: TImage;
    ImageHSt: TImage;
    procedure FormCreate(Sender: TObject);
    procedure savegraph1Click(Sender: TObject);
  private
    { Private declarations }
  public
    { Public declarations }
    deciaantal:integer;
    puntgrootte:integer;
    Ptmax:real;
    U14max:real;
    Ftmax:real;
    Gtmax:real;
    Htmax:real;
    Atmax:real;
    HStmax:real;
    Ptgraph:Tbitmap;
    U14graph:Tbitmap;
    Ftgraph:Tbitmap;
    Gtgraph:Tbitmap;
    Htgraph:Tbitmap;
    Atgraph:Tbitmap;
    HStgraph:Tbitmap;
    procedure TemporarySave;
    procedure InitiateGraphs;
    procedure DrawAxes(xxmin,xxmax,xxstap,yymin,yymax,yystap:real;ydeci:integer;var ggraph:Tbitmap;xxstring,yystring:string);
    procedure DrawPoint(var xn,yn,xnn,ynn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean;var ggraph:Tbitmap);
    procedure DrawpointPt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
    procedure DrawpointU14(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
    procedure DrawpointFt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
    procedure DrawpointGt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
    procedure DrawpointHt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
    procedure DrawpointAt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
    procedure DrawpointHSt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);

  end;

var
  Form1: TForm1;

implementation

{$R *.dfm}

uses unit2;

function known(r:real):boolean;
begin
  if ABS(r - unknown) < small
  then result := false
  else result := true;
end;


procedure Tform1.TemporarySave;
begin
  Ptgraph.savetofile(form2.logdir + 'Pt.bmp');           Ptgraph.free;  Ptgraph := nil;
  U14graph.savetofile(form2.logdir + 'U14.bmp');         U14graph.free; U14graph := nil;
  Ftgraph.SaveToFile(form2.logdir + 'Ft.bmp');           Ftgraph.Free;  Ftgraph := nil;
  Gtgraph.savetofile(form2.logdir + 'Gt.bmp');           Gtgraph.free;  Gtgraph := nil;
  Htgraph.savetofile(form2.logdir + 'Ht.bmp');           Htgraph.free; Htgraph := nil;
  Atgraph.savetofile(form2.logdir + 'At.bmp');           Atgraph.free; Atgraph := nil;
  HStgraph.savetofile(form2.logdir + 'HSt.bmp');         HStgraph.free; HStgraph := nil;
  form1.imagePt.picture.loadfromfile(form2.logdir + 'Pt.bmp');
  form1.imageU14.picture.loadfromfile(form2.logdir + 'U14.bmp');
  form1.imageFt.picture.loadfromfile(form2.logdir + 'Ft.bmp');
  form1.ImageGt.Picture.LoadFromFile(form2.logdir + 'Gt.bmp');
  form1.imageHt.picture.loadfromfile(form2.logdir + 'Ht.bmp');
  form1.imageAt.Picture.LoadFromFile(form2.logdir + 'At.bmp');
  form1.ImageHSt.Picture.LoadFromFile(form2.logdir + 'HSt.bmp');
end;


procedure TForm1.savegraph1Click(Sender: TObject);
var s:string;
    i:integer;
begin
  s := DateToStr(Date) + ' ' + Timetostr(Time);
  for i:= 1 to length(s) do if s[i] in ['/',':'] then s[i] := '_';
  case pagecontrol1.ActivePageindex of
   0 : imagePt.Picture.SaveToFile(form2.logdir + 'graph_'+ tabsheet1.caption + '_'+ s + '.bmp');
   1 : imageU14.Picture.SaveToFile(form2.logdir + 'graph_'+ tabsheet2.caption + '_'+ s + '.bmp');
   2 : imageFt.Picture.SaveToFile(form2.logdir + 'graph_'+ tabsheet3.caption + '_'+ s + '.bmp');
   3 : imageGt.Picture.SaveToFile(form2.logdir + 'graph_'+ tabsheet4.Caption + '_'+ s + '.bmp');
   4 : imageHt.Picture.SaveToFile(form2.logdir + 'graph_'+ tabsheet5.caption + '_'+ s + '.bmp');
   5 : imageAt.Picture.SaveToFile(form2.logdir + 'graph_'+ tabsheet6.Caption + '_'+ s + '.bmp');
   6 : imageHSt.Picture.Savetofile(form2.logdir + 'graph_'+ tabsheet7.caption + '_'+ s + '.bmp');
   end;
end;



procedure Tform1.InitiateGraphs;
var code :integer;
begin
  Ptmax := 10;
  U14max := 100;
  Ftmax := 1;
  Gtmax := 2;
  Htmax := 2;
  Atmax := 100;
  HStmax := 10;
  DrawAxes(0,maxduration,25,0,Ptmax,Ptmax/10,0,Ptgraph,'year','Population (x million)');
  DrawAxes(0,maxduration,25,0,U14max,U14max/10,0,U14graph,'year','Under 14 (in %)');
  DrawAxes(0,maxduration,25,0,Ftmax,Ftmax/10,1,Ftgraph,'year','Births (x million)');
  DrawAxes(0,maxduration,25,0,Gtmax,Gtmax/10,1,Gtgraph,'year','Growth rate (in %)');
  DrawAxes(0,maxduration,25,0,Htmax,Htmax/10,1,Htgraph,'year','Housewives (x million)');
  DrawAxes(0,maxduration,25,0,Atmax,Atmax/10,0,Atgraph,'year','Mean Age');
  DrawAxes(0,maxduration,25,0,HStmax,HStmax/10,0,HStgraph,'year','Household size');
end;




procedure TForm1.DrawAxes(xxmin,xxmax,xxstap,yymin,yymax,yystap:real;Ydeci:integer;var ggraph:Tbitmap;xxstring,yystring:string);
var xx,yy, xxas, yyas, i,l,l2,l3,imin,imax, hulpsize:integer;
    xi,yi:real;
    s:string;
begin
  ggraph := Tbitmap.create;
  hulpsize := 10;
  ggraph.canvas.Font.Size := hulpsize;
  ggraph.width := wwidth + 100;
  ggraph.height := hheight + 100;
  xxas := round(-xxmin/(xxmax-xxmin) * wwidth) + 50;
  if ((xxas < 0) or (xxas > ggraph.width)) then xxas := 0;
  yyas := round(yymax/(yymax-yymin)*hheight) + 50;
  if ((yyas < 0) or (yyas > ggraph.height)) then yyas := hheight;

  if form2.CheckBox1.Checked
  then ggraph.canvas.pen.color := clsilver
  else ggraph.Canvas.Pen.Color := clblack;

  //plaats horizontale lijnen
  imin := 0; while imin*yystap > yymin do dec(imin,1);
  imax := 0; while imax*yystap < yymax do inc(imax,1);
  for i:= imin +1 to imax+5
  do begin
     yi := i*yystap;
     xx := round(-xxmin/(xxmax-xxmin) * wwidth) + 50;
     yy := round(-(yi-yymax)/(yymax-yymin)*hheight) + 50;

     if form2.checkbox1.checked
     then begin
          ggraph.canvas.moveto(0,yy+1);
          ggraph.canvas.lineto(wwidth+100,yy+1);
          end
     else begin
          ggraph.Canvas.MoveTo(xxas,yy+1);
          ggraph.Canvas.lineto(xxas +10,yy+1);
          end;
     str(yi:3:ydeci,s);
     l := ggraph.Canvas.textwidth(s) + puntgrootte;
     if ((xx > 0) and (xx < ggraph.width))
     then ggraph.canvas.textout(xx - l,yy-7,s)
     else ggraph.canvas.textout(0 - l,yy-7,s);
     end;

  // plaats verticale lijnen
  imin := 0; while imin*xxstap > xxmin do dec(imin,1);
  imax := 0; while imax*xxstap < xxmax do inc(imax,1);
  for i:= imin to imax+10
  do begin
     xi := i*xxstap;
     xx := round((xi-xxmin)/(xxmax-xxmin) * wwidth) + 50;
     yy := round(yymax/(yymax-yymin)*hheight) + 50;

     if form2.checkbox1.checked
     then begin
          ggraph.canvas.moveto(xx,0);
          ggraph.canvas.lineto(xx,hheight+100);
          end
     else begin
          ggraph.Canvas.MoveTo(xx,yyas);
          ggraph.Canvas.LineTo(xx,yyas-10);
          end;
     str(xi+ form2.tbegin:3:0,s);
     l2 := ggraph.canvas.textwidth(s);
     if ((yy > 0) and (yy < ggraph.height))
     then ggraph.canvas.textout(xx - l2 div 2,yy+puntgrootte,s)
     else ggraph.canvas.textout(xx - l2 div 2,hheight,s);
     end;

  // plaats de twee assen
  ggraph.canvas.pen.color := clblack;
  ggraph.canvas.moveto(xxas,0);
  ggraph.canvas.lineto(xxas,yyas+1 {hheight+100});
  ggraph.canvas.moveto(0,yyas+1);
  ggraph.canvas.lineto(wwidth+100,yyas+1);

  //plaats de grootheden bij de assen
  ggraph.canvas.font.size := hulpsize + 2;
  {ggraph.Canvas.Font.Style := [fsbold];}
  ggraph.Canvas.Font.Orientation := 900;
  l3 := ggraph.Canvas.TextHeight(yystring);
  ggraph.canvas.textout(xxas -l - l3- puntgrootte,hheight + 30 {* 14 div 15},yystring);
  ggraph.Canvas.Font.Orientation := 0;
  ggraph.canvas.TextOut(wwidth div 4,yyas + puntgrootte + l3,xxstring);
  ggraph.canvas.font.size := hulpsize;
  {ggraph.Canvas.Font.Style := [];}
end;


procedure Tform1.DrawPoint(var xn,yn,xnn,ynn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean;var ggraph:Tbitmap);
begin
  if ((xnn > 0) and (xnn < ggraph.width) and (ynn > 0) and (ynn < ggraph.height))
  then begin
       {grafiek.canvas.pixels[xn,yn] := kleur; }
       ggraph.canvas.brush.color := kleur;
       ggraph.canvas.brush.style := bssolid;
       ggraph.canvas.pen.color := kleur;
       If ((verbinding) and (xn > 0) and (xn < ggraph.width) and (yn > 0) and (yn < ggraph.height))
       then begin
            ggraph.canvas.moveto(xn,yn);
            ggraph.canvas.lineto(xnn,ynn);
            end;
       ggraph.canvas.ellipse(xnn-grootte*puntgrootte,ynn-grootte*puntgrootte,xnn+grootte*puntgrootte,ynn+grootte*puntgrootte);
       xn := xnn;
       yn := ynn;
       end;
end;


procedure Tform1.DrawpointPt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
var Pnn,tnn:integer;
begin
  if known(tt) then tnn := round((tt- form2.tbegin)/maxduration * wwidth) + 50 else tnn := unknown;
  if known(PP) then Pnn := round(-(PP-Ptmax)/Ptmax * hheight) + 50 else Pnn := unknown;
  DrawPoint(tn,Pn,tnn,Pnn,kleur,grootte,verbinding,Ptgraph);
end;

procedure Tform1.DrawpointU14(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
var Pnn,tnn:integer;
begin
  if known(tt) then tnn := round((tt- form2.tbegin)/maxduration * wwidth) + 50 else tnn := unknown;
  if known(PP) then Pnn := round(-(PP-U14max)/U14max * hheight) + 50 else Pnn := unknown;
  DrawPoint(tn,Pn,tnn,Pnn,kleur,grootte,verbinding,U14graph);
end;

procedure Tform1.DrawpointFt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
var Pnn,tnn:integer;
begin
  if known(tt) then tnn := round((tt- form2.tbegin)/maxduration * wwidth) + 50 else tnn := unknown;
  if known(PP) then Pnn := round(-(PP-Ftmax)/Ftmax * hheight) + 50 else Pnn := unknown;
  DrawPoint(tn,Pn,tnn,Pnn,kleur,grootte,verbinding,Ftgraph);
end;

procedure Tform1.DrawpointGt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
var Pnn,tnn:integer;
begin
  if known(tt) then tnn := round((tt- form2.tbegin)/maxduration * wwidth) + 50 else tnn := unknown;
  if known(PP) then Pnn := round(-(PP-Gtmax)/Gtmax * hheight) + 50 else Pnn := unknown;
  DrawPoint(tn,Pn,tnn,Pnn,kleur,grootte,verbinding,Gtgraph);
end;

procedure Tform1.DrawpointHt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
var Pnn,tnn:integer;
begin
  if known(tt) then tnn := round((tt- form2.tbegin)/maxduration * wwidth) + 50 else tnn := unknown;
  if known(PP) then Pnn := round(-(PP-Htmax)/Htmax * hheight) + 50 else Pnn := unknown;
  DrawPoint(tn,Pn,tnn,Pnn,kleur,grootte,verbinding,Htgraph);
end;

procedure Tform1.DrawpointAt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
var Pnn,tnn:integer;
begin
  if known(tt) then tnn := round((tt- form2.tbegin)/maxduration * wwidth) + 50 else tnn := unknown;
  if known(PP) then Pnn := round(-(PP-Atmax)/Atmax * hheight) + 50 else Pnn := unknown;
  DrawPoint(tn,Pn,tnn,Pnn,kleur,grootte,verbinding,Atgraph);
end;

 procedure Tform1.DrawpointHSt(PP,tt:real;var Pn,tn:integer;kleur:Tcolor;grootte:integer;verbinding:boolean);
var Pnn,tnn:integer;
begin
  if known(tt) then tnn := round((tt- form2.tbegin)/maxduration * wwidth) + 50 else tnn := unknown;
  if known(PP) then Pnn := round(-(PP-HStmax)/HStmax * hheight) + 50 else Pnn := unknown;
  DrawPoint(tn,Pn,tnn,Pnn,kleur,grootte,verbinding,HStgraph);
end;



procedure TForm1.FormCreate(Sender: TObject);
begin
  deciaantal := 2;
  puntgrootte := 2;
end;

end.






unit NelderMead;

interface

const
maxparameters = 4;
dSimplex = 0.01;
Cpar : array[1..maxparameters] of real = (2.5,1.6,2,5.5);

type parameters = array[1..maxparameters] of real;


type TSimplex = class(Tobject)
       vertices: array[0..maxparameters] of parameters;
       values: array[0..maxparameters] of real;
       Margin:real;
       IHigh:integer;
       INearhigh:integer;
       ILow :integer;
       Centre: parameters;
       CentreValue:real;
       Reflection: parameters;
       ReflValue:real;
       Expansion:parameters;
       ExpaValue:real;
       Contraction:parameters;
       ContrValue:real;
       alpha:real;
       beta:real;
       gamma:real;
       delta:real;
       maxpar:integer;
       constructor New(a,b,c,d:real; npar:integer);
       procedure StartSimplex;
       procedure CalculateCentre;
       procedure CalculateReflection;
       procedure CalculateExpansion;
       procedure CalculateContraction;
       function CalculateValue(P:Parameters):real; virtual;
       procedure CalculateMargin;
       procedure DetermineHighAndLow;
       procedure Replace(n:integer;P:parameters; v:real);
       procedure Shrink;
       function Iteration:real;
       end;

implementation

uses unit2;


procedure TSimplex.Replace(n: Integer; P: parameters;v:real);
var i:integer;
begin
  for i:= 1 to maxpar do Vertices[n,i] := P[i];
  values[n] := v;
end;

constructor TSimplex.New(a,b,c,d: Real; npar:integer);
begin
  inherited create;
  alpha := a;
  beta := b;
  gamma := c;
  delta := d;
  maxpar := npar;
end;



procedure TSimplex.StartSimplex;
var i,j:integer;
    p, q:real;
begin
   p := dSimplex * (sqrt(maxpar+1) + maxpar - 1)/(maxpar * sqrt(2));
   q := dSimplex * (sqrt(maxpar+1)      -1)/(maxpar * sqrt(2));
   for j:= 1 to maxpar do Vertices[0,j] := Cpar[j] - (p + (maxpar-1) * q)/(maxpar+1);
   for i:= 1 to maxpar
   do for j:= 1 to maxpar
      do if j=i
         then Vertices[i,j] := Vertices[0,j] + p
         else Vertices[i,j] := Vertices[0,j] + q;
   for i:= 0 to maxpar do values[i] := CalculateValue(Vertices[i]);
   DetermineHighAndLow;
   CalculateCentre;
end;

procedure TSimplex.CalculateMargin;
var i,j:integer;
    r:real;
begin
  Margin := 0;
  for i:= 0 to maxpar
  do begin
     r := 0;
     for j:= 1 to maxpar do r := r + sqr(Vertices[i,j] - Centre[j]);
     if sqrt(r) > Margin then Margin := sqrt(r);
     end;
end;

procedure TSimplex.CalculateCentre;
var i,j:integer;
begin
  for i:= 1 to maxpar
  do begin
     Centre[i] := 0;
     for j:= 0 to maxpar do if j <> Ihigh then Centre[i] := Centre[i] + Vertices[j,i]/maxpar;
     end;
  CentreValue := CalculateValue(Centre);
  CalculateMargin;
end;

procedure TSimplex.CalculateReflection;
var i:integer;
begin
  for i:= 1 to maxpar
  do Reflection[i] := (1 + alpha) * Centre[i] - alpha * Vertices[Ihigh,i];
  ReflValue := CalculateValue(Reflection);
end;

procedure TSimplex.CalculateExpansion;
var i:integer;
begin
  for i:= 1 to maxpar
  do Expansion[i] := gamma * Reflection[i] + (1-gamma) * Centre[i];
  ExpaValue := CalculateValue(Expansion);
end;

procedure TSimplex.CalculateContraction;
var i:integer;
begin
  for i:= 1 to maxpar
  do Contraction[i] := beta * Vertices[IHigh,i] + (1-beta) * Centre[i];
  ContrValue := CalculateValue(Contraction);
end;

procedure TSimplex.Shrink;
var i,j:integer;
begin
  for i:= 0 to maxpar
  do if i <> ILow
     then begin
          for j:= 1 to maxpar
          do Vertices[i,j] := delta * Vertices[i,j] + (1-delta) * Vertices[ILow,j];
          values[i] := CalculateValue(Vertices[i]);
          end;
end;

procedure TSimplex.DetermineHighAndLow;
var r:real;
    i:integer;
begin
  r := values[0];
  ILow := 0;
  for i:= 1 to maxpar
  do if values[i] < r
     then begin
          Ilow := i;
          r := values[i];
          end;
  r := values[0];
  IHigh := 0;
  for i:= 1 to maxpar
  do if values[i] > r
     then begin
          IHigh := i;
          r := values[i];
          end;
  r := values[ILow];
  INearHigh := ILow;
  for i:= 0 to maxpar
  do if ((values[i] > r) and  (i<> IHigh))
     then begin
          INearHigh := i;
          r := values[i];
          end;
end;


function TSimplex.Iteration:real;
begin
  CalculateReflection;
  If ReflValue < values[Ilow]
  then begin
       CalculateExpansion;
       If ExpaValue < values[Ilow]
       then Replace(IHigh,Expansion, ExpaValue)
       else Replace(IHigh,Reflection, ReflValue);
       end
  else if NOT(ReflValue > values[INearHigh])
       then Replace(IHigh,Reflection, ReflValue)
       else begin
            if NOT(ReflValue > values[IHigh])
            then Replace(IHigh,Reflection, ReflValue);
            CalculateContraction;
            if Not(ContrValue > values[IHigh])
            then Replace(IHigh,Contraction, ContrValue)
            else Shrink;
            end;
  DetermineHighAndLow;
  CalculateCentre;
  result := CentreValue;
end;


function TSimplex.CalculateValue(P:Parameters):real;
var i:integer;
begin
  with form2
  do begin
     kappa := P[1]/1000;
     labda := P[2]/100;
     PBegin := P[3];
     PtotalJan := PBegin;
     calculategrowthrateHouseholds;
     Householdsize := 5.87 - 0.821 * growthrateHouseholds * 100;
     TestSimulationAdaptingBirths;
     error := fit;
     end;
  result := form2.error;
end;




end.



object Form2: TForm2
  Left = 0
  Top = 0
  Caption = 'form2'
  ClientHeight = 599
  ClientWidth = 1193
  Color = clBtnFace
  Font.Charset = DEFAULT_CHARSET
  Font.Color = clWindowText
  Font.Height = -11
  Font.Name = 'Tahoma'
  Font.Style = []
  OldCreateOrder = False
  OnCreate = FormCreate
  PixelsPerInch = 96
  TextHeight = 13
  object Panel1: TPanel
    Left = 8
    Top = 8
    Width = 313
    Height = 129
    TabOrder = 0
    object Label1: TLabel
      Left = 16
      Top = 16
      Width = 174
      Height = 19
      Caption = 'mortality parameters'
      Font.Charset = DEFAULT_CHARSET
      Font.Color = clWindowText
      Font.Height = -16
      Font.Name = 'Tahoma'
      Font.Style = [fsBold]
      ParentFont = False
    end
    object Label2: TLabel
      Left = 123
      Top = 47
      Width = 22
      Height = 13
      Caption = 'beta'
    end
    object Label3: TLabel
      Left = 209
      Top = 47
      Width = 26
      Height = 13
      Caption = 'theta'
    end
    object Label4: TLabel
      Left = 16
      Top = 88
      Width = 101
      Height = 13
      Caption = 'Black Death mortality'
    end
    object Label7: TLabel
      Left = 16
      Top = 47
      Width = 88
      Height = 13
      Caption = 'pre-crisis mortality'
    end
    object Label13: TLabel
      Left = 159
      Top = 88
      Width = 11
      Height = 13
      Caption = '%'
    end
    object Label14: TLabel
      Left = 203
      Top = 88
      Width = 29
      Height = 13
      Caption = 'decay'
    end
    object Edit1: TEdit
      Left = 150
      Top = 44
      Width = 42
      Height = 21
      TabOrder = 0
      Text = '0.064'
    end
    object Edit2: TEdit
      Left = 238
      Top = 44
      Width = 41
      Height = 21
      TabOrder = 1
      Text = '-6.3'
    end
    object Edit3: TEdit
      Left = 123
      Top = 85
      Width = 30
      Height = 21
      TabOrder = 2
      Text = '47'
    end
    object Edit6: TEdit
      Left = 238
      Top = 85
      Width = 41
      Height = 21
      TabOrder = 3
      Text = '0.043'
    end
  end
  object Panel2: TPanel
    Left = 8
    Top = 401
    Width = 313
    Height = 198
    TabOrder = 1
    object Label8: TLabel
      Left = 16
      Top = 23
      Width = 85
      Height = 19
      Caption = 'simulation'
      Font.Charset = DEFAULT_CHARSET
      Font.Color = clWindowText
      Font.Height = -16
      Font.Name = 'Tahoma'
      Font.Style = [fsBold]
      ParentFont = False
    end
    object Label9: TLabel
      Left = 16
      Top = 88
      Width = 69
      Height = 13
      Caption = 'starting year: '
    end
    object Label10: TLabel
      Left = 16
      Top = 131
      Width = 149
      Height = 13
      Caption = 'population size at starting year'
    end
    object Label11: TLabel
      Left = 16
      Top = 56
      Width = 47
      Height = 13
      Caption = 'duration: '
    end
    object Label12: TLabel
      Left = 204
      Top = 130
      Width = 28
      Height = 13
      Caption = 'million'
    end
    object Label21: TLabel
      Left = 116
      Top = 56
      Width = 27
      Height = 13
      Caption = 'years'
    end
    object Label6: TLabel
      Left = 16
      Top = 168
      Width = 204
      Height = 13
      Caption = 'pre-plague Household size at starting year'
    end
    object Label25: TLabel
      Left = 247
      Top = 128
      Width = 11
      Height = 13
      Caption = 'N:'
    end
    object Edit7: TEdit
      Left = 171
      Top = 126
      Width = 27
      Height = 21
      TabOrder = 0
      Text = '1.7'
    end
    object Button1: TButton
      Left = 226
      Top = 8
      Width = 75
      Height = 25
      Caption = 'OK'
      TabOrder = 1
      OnClick = Button1Click
    end
    object SpinEdit1: TSpinEdit
      Left = 91
      Top = 88
      Width = 51
      Height = 22
      MaxValue = 0
      MinValue = 0
      TabOrder = 2
      Value = 1086
    end
    object SpinEdit2: TSpinEdit
      Left = 69
      Top = 52
      Width = 41
      Height = 22
      MaxValue = 0
      MinValue = 0
      TabOrder = 3
      Value = 500
    end
    object Edit4: TEdit
      Left = 226
      Top = 164
      Width = 28
      Height = 21
      TabOrder = 4
      Text = '5'
    end
    object Button2: TButton
      Left = 113
      Top = 15
      Width = 107
      Height = 25
      Caption = 'Best fit Nelder Mead'
      TabOrder = 5
      OnClick = Button2Click
    end
    object Button3: TButton
      Left = 149
      Top = 46
      Width = 153
      Height = 25
      Caption = 'Best fit Population iterative'
      TabOrder = 6
      OnClick = Button3Click
    end
    object Button4: TButton
      Left = 157
      Top = 77
      Width = 144
      Height = 25
      Caption = 'Best fit Households iterative'
      TabOrder = 7
      OnClick = Button4Click
    end
    object Spinedit5: TSpinEdit
      Left = 264
      Top = 126
      Width = 39
      Height = 22
      MaxValue = 0
      MinValue = 0
      TabOrder = 8
      Value = 20
    end
  end
  object Panel3: TPanel
    Left = 8
    Top = 143
    Width = 313
    Height = 118
    TabOrder = 2
    object Label5: TLabel
      Left = 16
      Top = 16
      Width = 247
      Height = 19
      Caption = 'nuptiality/fertility parameters'
      Font.Charset = DEFAULT_CHARSET
      Font.Color = clWindowText
      Font.Height = -16
      Font.Name = 'Tahoma'
      Font.Style = [fsBold]
      ParentFont = False
    end
    object Label15: TLabel
      Left = 16
      Top = 83
      Width = 175
      Height = 13
      Caption = '(1-proportion widows)/Birth interval:'
    end
    object Label19: TLabel
      Left = 16
      Top = 53
      Width = 179
      Height = 13
      Caption = 'mean wedding age females (in years)'
    end
    object Edit5: TEdit
      Left = 197
      Top = 80
      Width = 29
      Height = 21
      TabOrder = 0
      Text = '0.33'
    end
    object SpinEdit4: TSpinEdit
      Left = 201
      Top = 50
      Width = 43
      Height = 22
      MaxValue = 0
      MinValue = 0
      TabOrder = 1
      Value = 24
    end
  end
  object Panel4: TPanel
    Left = 336
    Top = 8
    Width = 849
    Height = 513
    TabOrder = 3
    object Label16: TLabel
      Left = 16
      Top = 16
      Width = 54
      Height = 19
      Caption = 'output'
      Font.Charset = DEFAULT_CHARSET
      Font.Color = clWindowText
      Font.Height = -16
      Font.Name = 'Tahoma'
      Font.Style = [fsBold]
      ParentFont = False
    end
    object Label17: TLabel
      Left = 16
      Top = 279
      Width = 101
      Height = 13
      Caption = 'percentage spinsters'
    end
    object Label18: TLabel
      Left = 16
      Top = 52
      Width = 43
      Height = 13
      Caption = 'variables'
    end
    object Memo1: TMemo
      Left = 16
      Top = 71
      Width = 825
      Height = 182
      Lines.Strings = (
        'Memo1')
      ScrollBars = ssVertical
      TabOrder = 0
    end
    object Memo2: TMemo
      Left = 16
      Top = 298
      Width = 825
      Height = 180
      Lines.Strings = (
        'Memo2')
      ScrollBars = ssVertical
      TabOrder = 1
    end
    object CheckBox1: TCheckBox
      Left = 712
      Top = 20
      Width = 57
      Height = 17
      Caption = '  grid'
      TabOrder = 2
    end
  end
  object Panel5: TPanel
    Left = 8
    Top = 273
    Width = 313
    Height = 122
    TabOrder = 4
    object Label20: TLabel
      Left = 16
      Top = 8
      Width = 188
      Height = 19
      Caption = 'growthrate households'
      Font.Charset = DEFAULT_CHARSET
      Font.Color = clWindowText
      Font.Height = -16
      Font.Name = 'Tahoma'
      Font.Style = [fsBold]
      ParentFont = False
    end
    object Label22: TLabel
      Left = 24
      Top = 48
      Width = 33
      Height = 13
      Caption = 'kappa:'
    end
    object Label23: TLabel
      Left = 140
      Top = 48
      Width = 30
      Height = 13
      Caption = 'labda:'
    end
    object Label24: TLabel
      Left = 247
      Top = 48
      Width = 11
      Height = 13
      Caption = 'N:'
    end
    object Edit8: TEdit
      Left = 63
      Top = 46
      Width = 54
      Height = 21
      TabOrder = 0
      Text = '0.0025'
    end
    object Edit9: TEdit
      Left = 176
      Top = 46
      Width = 51
      Height = 21
      TabOrder = 1
      Text = '0.016'
    end
    object SpinEdit3: TSpinEdit
      Left = 264
      Top = 46
      Width = 37
      Height = 22
      MaxValue = 0
      MinValue = 0
      TabOrder = 2
      Value = 0
    end
  end
end






object Form1: TForm1
  Left = 0
  Top = 0
  Caption = 'Form1'
  ClientHeight = 576
  ClientWidth = 1050
  Color = clBtnFace
  Font.Charset = DEFAULT_CHARSET
  Font.Color = clWindowText
  Font.Height = -13
  Font.Name = 'Tahoma'
  Font.Style = []
  Menu = MainMenu1
  OldCreateOrder = False
  OnCreate = FormCreate
  PixelsPerInch = 96
  TextHeight = 16
  object Panel1: TPanel
    Left = 0
    Top = 0
    Width = 1050
    Height = 25
    Align = alTop
    TabOrder = 0
    ExplicitWidth = 773
  end
  object PageControl1: TPageControl
    Left = 0
    Top = 25
    Width = 1050
    Height = 551
    ActivePage = TabSheet7
    Align = alClient
    TabOrder = 1
    ExplicitWidth = 773
    object TabSheet1: TTabSheet
      Caption = 'Population'
      ExplicitWidth = 765
      object ImagePt: TImage
        Left = 0
        Top = 0
        Width = 1042
        Height = 520
        Align = alClient
        ExplicitLeft = 424
        ExplicitTop = 256
        ExplicitWidth = 105
        ExplicitHeight = 105
      end
    end
    object TabSheet2: TTabSheet
      Caption = 'Under 14'
      ImageIndex = 1
      ExplicitWidth = 765
      object ImageU14: TImage
        Left = 0
        Top = 0
        Width = 1042
        Height = 520
        Align = alClient
        ExplicitLeft = 384
        ExplicitTop = 208
        ExplicitWidth = 105
        ExplicitHeight = 105
      end
    end
    object TabSheet3: TTabSheet
      Caption = 'Fertility'
      ImageIndex = 2
      ExplicitWidth = 765
      object ImageFt: TImage
        Left = 0
        Top = 0
        Width = 1042
        Height = 520
        Align = alClient
        ExplicitLeft = 472
        ExplicitTop = 296
        ExplicitWidth = 105
        ExplicitHeight = 105
      end
    end
    object TabSheet4: TTabSheet
      Caption = 'Growth rate'
      ImageIndex = 3
      ExplicitWidth = 765
      object Imagegt: TImage
        Left = 0
        Top = 0
        Width = 1042
        Height = 520
        Align = alClient
        ExplicitLeft = 272
        ExplicitTop = 232
        ExplicitWidth = 105
        ExplicitHeight = 105
      end
    end
    object TabSheet5: TTabSheet
      Caption = 'Households'
      ImageIndex = 4
      ExplicitWidth = 765
      object ImageHt: TImage
        Left = 0
        Top = 0
        Width = 1042
        Height = 520
        Align = alClient
        ExplicitLeft = 384
        ExplicitTop = 264
        ExplicitWidth = 105
        ExplicitHeight = 105
      end
    end
    object TabSheet6: TTabSheet
      Caption = 'mean age of married women'
      ImageIndex = 5
      ExplicitWidth = 765
      object ImageAt: TImage
        Left = 0
        Top = 0
        Width = 1042
        Height = 520
        Align = alClient
        ExplicitLeft = 336
        ExplicitTop = 72
        ExplicitWidth = 105
        ExplicitHeight = 105
      end
    end
    object TabSheet7: TTabSheet
      Caption = 'Household size'
      ImageIndex = 6
      ExplicitWidth = 765
      object ImageHSt: TImage
        Left = 0
        Top = 0
        Width = 1042
        Height = 520
        Align = alClient
        ExplicitLeft = -184
        ExplicitTop = -21
        ExplicitWidth = 765
      end
    end
  end
  object MainMenu1: TMainMenu
    Left = 608
    object file1: TMenuItem
      Caption = 'file'
      object savegraph1: TMenuItem
        Caption = 'save graph'
        OnClick = savegraph1Click
      end
    end
  end
end




program Project2;

uses
  Vcl.Forms,
  Unit2 in 'Unit2.pas' {Form2},
  Unit1 in 'Unit1.pas' {Form1},
  NelderMead in 'NelderMead.pas';

{$R *.res}

begin
  Application.Initialize;
  Application.MainFormOnTaskbar := True;
  Application.CreateForm(TForm2, Form2);
  Application.CreateForm(TForm1, Form1);
  Application.Run;
end.