7+͈C:)=_gSIMP12 CRCSIMP12/TDOC SIMP12/TPASINITCMD PASD"ECHO PASfECHO COM|CLST2 PAS MINX PASPASS PAS FFTRUE PASMLTEST PASMM PAS MM SP PASSP  ML PAS! ML .RF PAS0RF BCRC Ver 5.0 CTL-S pauses, CTL-C aborts --> FILE: SIMP12/T.DOC CRC = 38 63 --> FILE: SIMP12/T.PAS CRC = 4E 63 --> FILE: INITCMD .PAS CRC = 76 CE --> FILE: ECHO .PAS CRC = 36 B2 --> FILE: ECHO .COM CRC = C4 CE --> FILE: LST2 .PAS CRC = 66 78 --> FILE: MINX .PAS CRC = 16 85 --> FILE: PASS .PAS CRC = AA CD --> FILE: FFTRUE .PAS CRC = 3A 01 --> FILE: MLTEST .PAS CRC = 7D D1 --> FILE: MM .PAS CRC = 0E D5 --> FILE: MM . CRC = 7C 34 --> FILE: SP .PAS CRC = F2 8E --> FILE: SP . CRC = 2D 34 --> FILE: ML .PAS CRC = B1 73 --> FILE: ML . CRC = 86 D6 --> FILE: RF .PAS CRC = D7 4F --> FILE: RF . CRC = 7F 12 DONE SIMP Version 1.2 TURBO Pascal (tm) May 16, 1984 D. M. Fritz-Rohner Post Office Box 9080 Akron, Ohio 44305 Thi versio i derive fro th articl 'Fittin Curve t Data b Marc S Cacec an Willia P Cacheri foun i BYT Magazine Volum 9 Numbe 5 May 198 an fro Versio 1. o th Pascal/ڠ (tm implementation Se Not fo lis o conversio effects Thi versio i intende t ru i CP/ (tm compatibl environment i conjunctio wit TURB Pasca. Thi versio furthe reorganize an extend th referenc code Se als SIMP11/Z.DO fro SIMP11/Z.LBR. Th structur o SIM 1. Pascal/ i largel preserve wit additio o boundar conditio testin an modularizatio ofcomman lin argumen processing. Almos al젠 non-academi extremizatio problem hav constraint i th for o auxiliar relation an boundar conditions Linkag t generi routin name ffTes i pro vide t retur boolea tru i th functio bein extremize i feasibl fo th specifie parameter an false otherwise I procedur SIMPLE wher reflection expansion o contractio ar applie i searchin fo bette extremum th teste poin mus als b feasibl a reporte b ffTest Necessarily thi featur require tha th spac containe b th initia SIMPLEؠ b feasible Se Example Multi-Linea Fi fo illustrativ applicatio o thi feature. TURB Pasca Versio 1. overwrite par o th CP/ comman tai buffe wit jum vector Onl th firs 3 byte o th comman tai ar availabl there Thi i generall enoug fo tw ful siz fil descriptor an separator bu allow n growth Th parse code i INITCMD.PA parse th reduce CP/ TURB comman tail Othe version ar availabl fo TURB tha includ explici I/ redirection ambiguou fil nam processing an pipe usin th CC comman lin image Lac o thes feature i no necessaril handica i numerica tool environment a oppose t softwar tool environments ECH illustrates usage of INITCMD routines. TURB i organize i wa tha doe no us separat compilatio an therefor doe no requir a explici linker Th routin SIMP12/T.PA i structure t us generi module wher selectio o include routine bind differen extremiza tio criteri an fittin function wit feasibilit tes an comman lin parse t for particula build Th default i th distribute versio o SIM 1. TURB Pasca are leas square minimization th referenc articl Michaelis-Mente function a unconditiona stu returnin tru fo ffTest an th reduce comman tai parser. TURB Pasca Versio 2. i reporte t provid segmentatio an overla facilitie whic woul offe another perhap mor convenient approac t progra structur. USAGE Ente TURBO cal u th (W)or fil SIMP.PAS ente (E)di an䠠 enabl堠 ($I)nclusio o th堠 appropriat堠 extremize function extremizatio criterion an feasiblilit test Exi edit selec (C)om-fil i th (O)ption menu an (C)ompile (Q)ui TURB an execut th generate build Suppos tha th extremize functio i th referenc articl Michaelis-Mente function th extremizatio criterio i leas squares an th feasibilit tes i alway true The th selectio preambl o SIMP.PA look lik th followin afte edit; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} {--$I b:? parameterized function } {$I b:mm.pas Michaelis-Menten, see ref. } {--$I b:sp.pas simple odd polynomial } {--$I b:ml.pas multi-linear fit } {--$I b:rf.pas 'root finder', see SIMP.DOC } {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} {--$I b:? extremization criterion } {$I b:lst2.pas normalized least squares } {--$I b:minx.pas minimax } {--$I b:pass.pas unity gain } {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} {--$I b:? feasibility test } {$I b:ffTrue.pas always returns true } {--$I b:mlTest.pas fease multi-linear fit } {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} {--$I b:? command parser } {$I b:initcmd.pas CP/M TURBO command tail parse} {--$I b:altcmd.pas CP/M CCP command line parse } {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} The command; A>simp mm take contro dat fro th fil M an execute outpu defaulte t th CP/ con device or; A>simp mm lst: with output to the CP/M lst: device, or; A>simp mm mmlst2.dat fo outpu t fil MMLST2.DA o th defaul drive Consol outpu ma b mad explici with; A>simp mm con: EXAMPLES Fo thi distributio th basi modul i SIMP Th primitive ar QUIԠ an INITCM whic contain InitAr an GetArg Th extremizatio criteri ar LST2 MINX an PASS Th extremize function ar MM SP M an R wit correspondin contro dat files Th feasibilit routine ar FFTRUŠ an MLTEST. SIMP.PA contain th mai progra an th SIMPLE modules QUIT.SR contain jum t CP/M' war boo rela instruction. I additio t th leas square extremizatio criterion LST2.PAS thi distributio include minima routine MINX.PA t facilitat searc fo parameter tha produc fi suc tha th maximu absolut valu o th differenc betwee measure an estimate value i minimized Minima als present usefu extremu criterio wher differentiabilit i los an linea method ar generall no helpful. Th extremizatio criterio PASS.PA i intende fo appli cation wher a extremu i sough an th extremizatio cri terio itsel PASSe th funtio valu withou referenc t 'fitted dat set Se Example: 'Roo Finder below. EXAMPLES Th example provide wit thi distributio ar Michaelis- Menten simpl polynomial multi-linea fit an 'roo finder'. Example: Michaelis-Menten The file MM contains the data; 200 - maximum number of iterations 5e-3 - extremum criterion convergence 2 - number of parameters 0.2 1.0 1.0 - 'a' initial, step, and convergence 3.0 1.0 1.0 - 'b' " " " " 6 - number of data values 0.172 1.68 - dependent variable, independent 1 0.250 3.33 - " " " 2 0.286 5.00 - ... 0.303 6.67 - ... 0.334 10.0 - 0.384 20.0 - " " " 6 Notic tha thes value ar arrange differentl an th ste size an convergenc criteri ar differen fro th refer enc article. The configuration selected above, namely; {$I b:mm.pas Michaelis-Menten, see ref. } {$I b:lst2.pas normalized least squares } {$I b:ffTrue.pas always returns true } {$I b:altcmd.pas CP/M CCP command line parse } compose progra tha seek bes fi base o leas square fo th referenc function. The command; A>simp mm test put th result i fil TEST Thi fil contain 20 lines th las o whic is; 200 4.238095E-01 2.451936E+00 5.284773E-03 5.284773E-03 where 20 mean tha th convergenc criteri wer no satis fied th firs an secon number ar th 'a an 'b coeffi cient a th verte wit leas extremu valu give b th thir number an th las numbe i th averag SIMPLE extremu fo al vertices Thi say that o th average th fi i goo t abou 0.00 fo eac coordinate. The following configuration; {$I b:mm.pas Michaelis-Menten, see ref. } {$I b:minx.pas normalized least squares } {$I b:ffTrue.pas always returns true } {$I b:altcmd.pas CP/M CCP command line parse } result i progra tha seek 'bes fit base o minima criterio fo th referenc article' functio MM Th followin command; A>simp mm test performs the extremization and puts the results in the file TEST on the default drive. The last line of that file is; 200 4.241221E-01 2.470153E+00 6.501847E-03 6.501847E-03 where 200 is the maximum iteration number - meaning that the iteration did not converge - and the first two numbers are the coefficient a before an th las tw number ar th extremu functio a th verte wit th bes fi an th averag fo al vertices respectively Th bes fi show minima erro 0.00. Selectio o extremu functio demand som insigh t th problem Fo th referenc proble i depend o th variance o th experimen an measuremen processe i linea estima tio sense Straigh forwar computationa procedure lik leas square o minima ma b equall wrong. Suppos yo hav choic betwee tw simila representa tion an kno tha on ha a averag erro abou equa t th other' maximu error Whic to choose? Example: Simple Polynomial Th pape ' Simpl Minima Algorithm b Steve A Ruzinsk foun i Dr Dobb' Journal Numbe 93 July 198 illustrate th propose algorith wit MBASI (tm compatibl progra tha use thre ter od polynomia t approximat th SIΠ functio ove th interva t abou pi/2 Se FLUF11.DOC. The data file SP contains corresponding information; 100 1.e-4 3 1.0 1.0 1.0 -.167 1.0 1.0 .00833 1.0 1.0 50 3.141076E-02 2.000000E-02 6.279052E-02 4.000000E-02 : : 9.995066E-01 9.800000E-01 1.000000E+00 1.000000E+00 wher th coefficient ar initialize t th approximat value for the corresponding Taylor series. Th followin configuratio selection produc leas square approximation; {$I b:sp.pas simple odd polynomial } {$I b:lst2.pas normalized least squares } {$I b:ffTrue.pas always returns true } {$I b:altcmd.pas CP/M CCP command line parse } which, when compiled and executed with the command; A>simp sp test take 14 second t produc th fil TES wit on hundre an one lines, the last of which is; 100 1.571301E+00 -7.02287E-01 1.380962E-01 5.112938E-03 6.206377E-03 Th followin include for configuratio t comput minima fi i thi problem; {$I b:sp.pas simple odd polynomial } {$I b:minx.pas minimax } {$I b:ffTrue.pas always returns true } {$I b:altcmd.pas CP/M CCP command line parse } Th command; A>simp sp test execute i 12 second an produce th dat fil TES wit th las line; 100 1.548642E+00 -6.22970E-01 8.073769E-02 8.232326E-03 8.428058E-03 whic show minima erro o abou .008 Suppos w us thi estimat t restar SIMP Wha happens Th followin show th firs an las line o th resultin dat fil produce i 10 seconds; 0 1.548642E+00 -6.22970E-01 8.073769E-02 8.232576E-03 7.568654E-01 : : 94 1.570308E+00 -6.41893E-01 7.167449E-02 9.655501E-05 9.721492E-05 whic show tha th specifie 1.e- convergenc toleranc wa achieved Th value reporte fo FLU ar roughly take fro Figur o th article; 1.5706264 -0.64322566 0.07270740 sho reasonabl agreemen considerin tha n qualificatio o anybody' SIΠ functions floatin poin conversion o compu tation i established. SIMP' slo convergenc i produce b th degre o con tractio require t maneuve i th minima space' cus geo metry Restartin SIM allow ne portio o th extremu spac t b explore befor th SIMPLE figur contracts Ther als come poin whe th SIMPLE i s diminishe tha th finit mathematic o th machin representatio o th extremu spac produc degenerat simple whe som dimensio vanishe o th finit numbe o accessibl point ar al unacceptable. SIM need anothe trick Th proble i ho t projec th SIMPLEؠ hyperhedro int reduce spac wit it geometr conforma t th locall reduce space. Example: Milti-Linear Fit The data in file ML are; 200 - maximum number of iterations 1.e-4 - extremum criterion convergence 6 - number of parameters 0.00000 0.1 1.0 - breakpoint coordinates, initial 0.00000 0.1 1.0 - steps and convergence criteria 0.30000 0.1 1.0 0.29552 0.1 1.0 0.60000 0.1 1.0 0.56464 0.1 1.0 7 - number of data values 0.00000 0. - dependentindependen variable 1 0.09983 0. - 2 0.198669 0.2 0.295520 0.3 0.389418 0.4 0.479426 0.5 0.564642 0.6 The configuration selected is shown by the following; {$I b:ml.pas multi-linear approximatin } {$I b:minx.pas minimax } {$I b:mlTest.pas multi-linear feasibility test} {$I b:altcmd.pas CP/M CCP command line parse } compose progra tha seek bes fi base o minima fo multi-linea approximatio fo SI tabl o th interva t abou pi/6 Th approximatio comprise tw lin segment wit initia breakpoint se a th first middle an las dat points respectively Th feasibilit tes doe no requir th breakpoint t sta i th range define b th dat points onl tha th independen variabl value remai i ascendin order Th breakpoin functio evaluatio set th functio valu t th respectiv endpoin valu i th independen vari abl i outsid th breakpoin range. Th command; A>simp ml test execute i 8 second an put th result i fil TEST Thi fil comprise 20 lines th first and last of which are; 0 0.000000E+00 0.000000E+00 3.000000E-01 2.955200E-01 6.000000E-01 5.646400E-01 4.492667E-03 8.078624E-02 : : 200 2.012359E-03 1.677780E-03 2.940030E-01 2.917680E-01 5.972280E-01 5.663182E-01 1.677916E-03 1.678885E-03 Therefore, the breakpoints are roughly; (.0020124,.0016778) (.2940030,.2917680) (.5972280,.5663182) with a minimax fit of about .00168. Example: 'Root Finder' In this example the routine RF.PAS codes a form to illus trate SIMPLEX application to extremization problems without fit conditions, functions that are not differentiable, and have a cusp geometry and multiple extrema. Suppose a second order polynomial with real roots at (x = -1.0) and (x = 2.0); 2 f(x) = (x + 1) * (x - 2) = x - x - 1 I thi exampl th parameter ar se t approximat value o th root an th independen dat value ar se t th polynomia coefficients B definition whe th valu o get t roo th polynomia evaluate t zero I th absolut valu o th polynomia i taken the th zer crossing becom minima Not tha thes minim ar cusp wit discontinuou firs derivatives. Further sinc multipl root ar bein sought th su o th absolut value o th polynomia evaluate a eac roo estimateca b mad th extremize functio returne t th SIMPLE procedure through PASS.PAS; m _____ \ ff = \ abs [ (..(x * a + x ) * a + ..) + x ] / 1 i 2 i n /____ i = 1 where; a[1..m] become the roots, x[1..n] become the coefficients, and y[1..n] dummied for consistency Th fil R contain th data; 100 - number of iterations 1e-5 - extremum convergence criterion 2 - number of 'parameters' 0.0 1.0 1e-5 - initial 'parameter' values, 0.0 1.0 1e-5 - steps, and convergence values 3 - number of 'data points' 0.0 1.0 - data point values 0.0 -1.0 - 0.0 -2.0 - The configuration produced by the include selections; {$I b:rf.pas 'root finder', see SIMP.DOC } {$I b:pass.pas unity gain } {$I b:ffTrue.pas always returns true } {$I b:altcmd.pas CP/M CCP command line parse } and performed with the command; A>simp rf test executes a program that seeks a best fit based on the specified extremu an put th result i fil TEST Th comman execute in 19 seconds and the file contains 80 lines, the last of which is; 79 -1.00000E+00 -1.00000E+00 5.671813E-06 9.414271E-06 where; 79 is the iteration number at which the convergence cri teria were jointly satisfied. Note that this function has extrema at (-1.0,-1.0) and (2.0,2.0) in addition to the 'correct' minimum at (-1.0,2.0). Dependin o initia conditions SIM cheerfull scramble int th neares cus valle an proceed downhil t th neares extremum. This example illustrates certain features of a class of extremization problems. It is not proposed for finding roots. Its defects also include problems with complex roots, osculating saddle points, and appropriate definition of 'abs' for complex values in addition to multiple extrema and slow convergence. There is no substitute for judgement or experience. Each algorithmic conclusion must be tested for reasonableness. McComb's Only Commandment states: 'Arrange the programming so that intelligence can successfully interfere with the computing.' NOTES 1. Conversion showed the following distinctions; i. String declaration syntax slightly different, e.g., var str : string 14 ; { Pascal/Z } var str : string[14] ; { TURBO } ii. Constant folding slightly different, e.g., const mxp1 = mx + 1 ; { Pascal/Z } const mxp1 = 9 ; { TURBO } iii. TURBO more insistent on Pascal scope delimitation, e.g.; if ... then begin ... end ; else begin ... end ; is acceptable to Pascal/Z but is rejected by TURBO unless the semicolon before the else is removed. iv. Association of file names and file variables is slightly different, e.g., reset(fvar,fname) ; { Pascal/Z } assign(fvar,fname) ; { TURBO } reset(fvar) ; v. String building is slightly different; append(sname,' ') ; { Pascal/Z } sname := concat(sname,' ') ; { TURBO } vi. TURBO provides explicit extensions to access environment. To get a high level exit and provoke warm boot, QUIT was programmed as an external procedure for Pascal/Z but was coded with a BIOS procedure call in TURBO. v. File open checks slightly different; eof(...) { Pascal/Z } IOresult { TURBO } 2. An interesting situation arises when comparing TURBO Pascal performance with compile-assemble-link translators like Pascal/Z. Suppose that comparable modular structures have been developed. The development cycle comparison is between TURBO's edit-(save)-compile-(load)-go and Pascal/Z's edit- compile-assemble-link-go. However, the 'correct' runtime comparison is between TURBO's compile and go and Pascal/Z's link, load and go. CP/M (tm) Digital Research Pascal/Z (tm) Ithaca Intersystems, (c) Jeff Moskow TURBO (tm) Borland International MBASIC (tm) MicroSoft UNIX (tm) AT&T Information Systems  { Copyright (C) 1984 D.M. Fritz-Rohner } {$I- suppress default I/O error } program main ; { } { Name: SIMPLEX - Apply SIMPLEX algorithm. } { } { Version: 1.2/Turbo (tm) Date: 1984 July 7 } { } { Purpose: Search for parameters to extremize } { specified criterion using specified function. } { } { This program is a generic application module } { intended to be particularized by selection or } { insertion of includes to bind routines for } { specified extremization criteria and specified } { functions. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { Reference: 'Fitting Curves to Data' } { Marco S. Caceci and William P. Cacheris } { BYTE Magazine, Volume 9, Number 5 } { May, 1984, pp. 340-362 } { } const mx = 8 ; { Max number of parameters } mxp1 = 9 ; { Number of SIMPLEX vertices } nx = 64 ; { Maximum number of data } ncxx = 1024 ; { Max number of itertions } CONTRACT = 0.618 ; { Contraction coefficient } REFLECT = -1.0 ; { Reflection coefficient } EXPAND = 1.618 ; { Expansion coefficient } BLANK = ' ' ; MAXLINE = 31 ; MAXSTR = 14 ; MAXARG = 8 ; type vertex = array[1..mx] of real ; vector = array[1..nx] of real ; index = 0..255 ; argstr = string[MAXSTR] ; argvec = array[1..MAXARG] of argstr ; var Narg : integer ; Argv : argvec ; fin : text ; Infile : string[14] ; fout : text ; { Input,output file vars } Outfile : string[14] ; err : text ; { Error file var } m : integer ; { Number of parameters } mi : integer ; mp1 : integer ; { m + 1 } simp : array[1..mxp1] of vertex ; { Vertices } Ssimp : array[1..mxp1] of real ; { Extrema } csimp : vertex ; { Centroid } dsimp : vertex ; { Parameter step sizes } asimp : vertex ; { Avg SIMPLEX vertex values } aSsimp : real ; { Avg SIMPLEX extremum value } cca : vertex ; { Parameter convergence } tsimp : vertex ; { Test vertex } Stsimp : real ; { Extremum value at test vertex } ccS : real ; { Extremum convergence } n : integer ; { Number of data } ni : integer ; ys : vector ; { Dependent data vector } xs : vector ; { Independent data vector } nc : integer ; { Iteration count } ncx : integer ; { Max number of iterations } i,j : index ; ix : index ; inn : index ; { Index to vertex with least } inx : index ; { Index to vertex with greatest } done : boolean ; { Termination flag } a : real ; { Intermediate variable } c : real ; NULL : real ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure quit ; begin bios(0) { warm boot } end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} {--$I b:? parameterized function } {$I b:mm.pas Michaelis-Menten, see ref. } {--$I b:sp.pas simple odd polynomial } {--$I b:ml.pas multi-linear approximation } {--$I b:rf.pas 'root finder', see SIMP.DOC } {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} {--$I b:? extremization criterion } {$I b:lst2.pas normalized least squares } {--$I b:minx.pas minimax } {--$I b:pass.pas unity gain } {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} {--$I b:? feasibility test } {$I b:ffTrue.pas always returns true } {--$I b:mlTest.pas test multi-linear fit } {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} {$I b:initcmd.pas CP/M TURBO command tail parse} {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure HuhEOF ; begin writeln(err,'Premature EOF.') ; quit end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure SimpData ; { } { Reads data from disk file FILENAME. Terminate file } { with EOF immediately after last digit. Data in the } { order; } { } { ncx - maximum number of iterations } { ccS - convergence criterion } { m - number of parameters } { a(1) da(1) cca(1) - initial parameter values, } { a(2) da(2) cca(2) steps sizes, and convergence } { a(3) da(3) cca(3) criteria, } { : : : } { a(m) da(m) cca(m) } { n - number of points } { ys(1) xs(1) - data values } { ys(2) xs(2) } { ys(3) xs(3) } { : : } { ys(n) xs(n) } { } begin read(fin,ncx) ; { Read maximum iteration count } if eof(fin) then HuhEOF ; if ncx < 1 then begin writeln(err,'Unreasonable Maximum Iteration Count: < 1') ; quit end ; if ncx > ncxx then begin writeln(err,'Unexpected Maximum Iteration Count: > ',ncxx:4) ; ncx := ncxx end ; read(fin,ccS) ; { Read convergence criterion } if eof(fin) then HuhEOF ; if ccS <= 0.0 then begin writeln(err,'Unreasonable Convergence Criterion: <= 0.') ; quit end ; read(fin,mi) ; { Read number of parameters } if eof(fin) then HuhEOF ; if mi < 1 then begin writeln(err,'Unreasonable Parameter Order: < 1') ; quit end ; if mi > mx then begin write(err,'Unexpected Parameter Order: m > ',mx:2) ; writeln(err,' : Excess Ignored') ; m := mx ; end else m := mi ; for j := 1 to m do begin { Read initial parameter values } read(fin,simp[1,j],dsimp[j],cca[j]) ; if eof(fin) then HuhEOF ; if dsimp[j] = 0.0 then begin writeln(err,'Unreasonable Parameter Step Size: = 0.0') ; quit end ; if cca[j] <= 0.0 then begin writeln(err,'Unreasonable Convergence Parameter: <= 0.') ; quit end end ; if mi > mx then for j := mx+1 to mi do begin { Gobble excess parameter values } read(fin,NULL,NULL,NULL) ; if eof(fin) then HuhEOF end ; read(fin,ni) ; { Read number of points } if ni > 0 then begin if eof(fin) then HuhEOF ; if ni < m then begin writeln(err,'Unreasonable Data Order: n < m') ; quit end ; if ni > nx then begin write(err,'Unexpected Data Order: n > ',nx:3) ; writeln(err,' : Excess Ignored') ; n := nx ; end else n := ni ; for i := 1 to n do begin { Read initial parameter values } read(fin,ys[i],xs[i]) ; if ((i < n) and eof(fin)) then HuhEOF end ; if ni > nx then for i := nx+1 to ni do begin { Gobble excess data values } read(fin,NULL,NULL) ; if ((i < ni) and eof(fin)) then HuhEOF end end else n := 0 end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure IndExt ; { } { Find indices for vertices with least and greatest } { values of extremum function. } { } begin inn := 1 ; inx := 1 ; for i := 2 to mp1 do if Ssimp[i] < Ssimp[inn] then inn := i else if Ssimp[i] > Ssimp[inx] then inx := i end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure SimpInit ; { } { Initialize SIMPLEX vertices, extremum vector, } { averages, and extrema indices. } { } begin mp1 := m + 1 ; { Compute number of SIMPLEX vertices } Ssimp[1] := SS(simp[1]) ; { Compute extremum value } aSsimp := Ssimp[1] ; for i := 2 to mp1 do begin { Initialize simplex vertices } for j := 1 to m do if j = i - 1 then simp[i,j] := simp[1,j] + dsimp[j] else simp[i,j] := simp[1,j] ; Ssimp[i] := SS(simp[i]) ; { Initialize extremum values } aSsimp := aSsimp + Ssimp[i] end ; aSsimp := aSsimp / mp1 ; { Initialize averages } for j := 1 to m do begin a := 0.0 ; for i := 1 to mp1 do a := a + simp[i,j] ; asimp[j] := a / mp1 end ; IndExt end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure Center ; { } { Compute center of SIMPLEX hyperplane face not con- } { taining vertex with greatest extremum value. } { } begin for j := 1 to m do begin c := 0.0 ; for i := 1 to mp1 do if i <> inx then c := c + simp[i,j] ; csimp[j] := c / m end end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure Offset ; { } { Compute difference vector from centroid to vertex } { with greatest extremum value. } { } begin for j := 1 to m do dsimp[j] := simp[inx,j] - csimp[j] end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure Project ( pc : real ) ; { } { Compute test vertex using centroid, Offset, and pro-} { jection coefficient and compute value of extremum } { function. } { } begin for j := 1 to m do tsimp[j] := csimp[j] + pc * dsimp[j] ; Stsimp := SS(tsimp) end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure Replace ( it : index ) ; { } { Substitute test vertex for specified vertex. } { } begin for j := 1 to m do simp[it,j] := tsimp[j] ; Ssimp[it] := Stsimp end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure SHRINK ; { } { Shrink SIMPLEX with respect to vertex with least } { value of extremum function. } { } begin for j := 1 to m do begin a := simp[inn,j] ; for i := 1 to mp1 do if i <> inn then simp[i,j] := a + CONTRACT * (simp[i,j] - a) end end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure Simplex ; { } { Perform one iteration using SIMPLEX projections to } { explore parameter space seeking extremum. } { } begin IndExt ; { Find indices for least and greatest } Center ; { Compute centroid vector } Offset ; { Compute difference vector } Project (REFLECT) ; { Try reflection } if ffTest(tsimp) and (Stsimp <= Ssimp[inn]) then begin Replace (inn) ; { Accept reflection and .. } Project (EXPAND) ; { Try expansion } if ffTest(tsimp) and (Stsimp <= Ssimp[inn]) then Replace (inn) { Accept expansion } end else begin { Try Replacement } if ffTest(tsimp) and (Stsimp <= Ssimp[inx]) then Replace (inx) { Accept Replacement } else begin Project (CONTRACT) ; { Try contraction } if ffTest(tsimp) and (Stsimp <= Ssimp[inx]) then Replace (inx) else SHRINK { Try shrinkage } end end end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure Results ; begin write(fout,nc:4) ; { Show Results } for j := 1 to m do write(fout,' ',simp[inn,j]:12) ; writeln(fout,' ',Ssimp[inn]:12,' ',aSsimp:12) end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} begin { Main program } nc := 0 ; arginit (Narg,Argv) ; { Parse command line and open files } assign(err,'con:') ; rewrite(err) ; { open error message file } if narg < 1 then begin writeln(err,'Usage: simp Infile [Outfile]') ; quit end ; if getarg(1,Infile,MAXARG) then begin assign(fin,Infile) ; reset(fin) ; { open Infile for input } if IOresult <> 0 then begin writeln(err,'Unable to Open ',Infile,' for Input.') ; quit end end ; if getarg(2,Outfile,MAXARG) then begin assign(fout,Outfile) ; rewrite(fout) ; { open Outfile for output } if IOresult <> 0 then begin writeln(err,'Unable to Open ',Outfile,' for Output.') ; quit end end else begin assign(fout,'con:') ; rewrite(fout) { open console for default output } end ; SimpData ; { Read control and data vectors } SimpInit ; { Initialize SIMPLEX } Results ; { Show initial parameters } repeat { Main loop } Simplex ; { Perform one iteration } done := true ; for j := 1 to m do begin { Use Cauchy parameter convergence test } a := 0.0 ; for i := 1 to mp1 do a := a + simp[i,j] ; a := a / mp1 ; if abs(asimp[j] - a) > cca[j] then done := false ; asimp[j] := a end ; aSsimp := 0.0 ; { Use absolute extremum convergence test } for i := 1 to mp1 do aSsimp := aSsimp + Ssimp[i] ; aSsimp := aSsimp / mp1 ; if abs(aSsimp) > ccS then done := false ; nc := nc + 1 ; { Increment iteration count } Results { Show iteration Results } until (done or (nc = ncx)) ; close(fout) end.  { Copyright (C) 1984 D.M. Fritz-Rohner } procedure ArgInit ( var narg : integer; var argv : argvec ) ; { } { Name: ArgInit - ARGument vector INITialization. } { } { Version: 1.0/TURBO (tm) Date: 1984 May 16 } { } { Purpose: Parse reduced command tail just enough to } { extract argument strings. Note that argument } { counting conforms to Pascal's one origin indexing } { so that argument indices run from 1 to narg. } { } { Notice that command tail length is limited because } { TURBO writes a jump vector over the command buffer } { beginning at offset 20H. Therefore, only 31 bytes } { of the command tail are available in default line } { buffer. See ALTCMD.PAS if the entire command line } { must be parsed. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } const MAXLINE = 31 ; type xstrptr = ^xstr ; xstr = record len : char ; tail : array[1..MAXLINE] of char ; end ; var arg : argstr ; { intermediate variable } cptr : xstrptr ; { command tail pointer } icl : integer ; { command tail index } ncl : integer ; { number command chars } begin cptr := ptr(128) ; { command line tail } ncl := ord(cptr^.len) ; if (ncl > MAXLINE) then ncl := MAXLINE ; narg := 0 ; icl := 1 ; while ((icl <= ncl) and (narg < MAXARG)) do begin arg := '' ; { gobble leading blanks } while ((icl <= ncl) and (cptr^.tail[icl] = BLANK)) do icl := icl + 1 ; { extract argument } while ((icl <= ncl) and (cptr^.tail[icl] <> BLANK) and (length(arg) < MAXSTR)) do begin arg := concat(arg,cptr^.tail[icl]) ; icl := icl + 1 end ; { skip argument tail if truncated } while ((icl <= ncl) and (cptr^.tail[icl] <> BLANK)) do icl := icl + 1 ; { update argument vector } if (length(arg) > 0) then begin narg := narg + 1 ; argv[narg] := arg end end end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} function GetArg (n : integer; var s : argstr; maxstr : integer) : boolean ; { } { Name: GetArg - GET and copy command line ARGument. } { } { Version: 1.0/Turbo (tm) Date: 1984 May 16 } { } { Purpose: Copy n-th argument from command line argu- } { ment vector to string s. Note that maxstr is dummy. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { Reference: 'Software Tools in Pascal' } { Brian W. Kernighan and P. J. Plauger } { Addison-Wesley Publishing Company, 1981 } { pp. 350 } { } begin if (n > 0) and (n <= Narg) then begin s := Argv[n] ; GetArg := true end else GetArg := false end ;  { Copyright (C) 1984 D.M. Fritz-Rohner } program main ; { } { Name: ECHO - ECHO command line arguments. } { } { Version: 1.0/TURBO (tm) Date: 1984 May 16 } { } { Purpose: Demonstrate command tail parsing using } { CP/M (tm) line buffer contents. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } const BLANK = ' ' ; MAXLINE = 31 ; MAXSTR = 14 ; MAXARG = 8 ; type argstr = string[MAXSTR] ; argvec = array[1..MAXARG] of argstr ; var Narg : integer ; Argv : argvec ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} {$I b:initcmd.pas parse command line } {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure echo ; { } { Name: ECHO - ECHO command line arguments. } { } { Version: 1.0/Turbo (tm) Date: 1984 May 16 } { } { Purpose: Echo command line arguments to output. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { Reference: 'Software Tools in Pascal' } { Brian W. Kernighan and P. J. Plauger } { Addison-Wesley Publishing Company, 1981 } { pp. 45 et seq } { } var i : integer ; arg : argstr ; begin i := 1 ; while (GetArg(i,arg,MAXSTR)) do begin writeln(output,arg) ; i := i + 1 end end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} begin { main program } ArgInit (Narg,Argv) ; { form argument list } echo { display argument list } end. ͫCopyright (C) 1984 BORLAND IncA SDS VDB-8024selectedP=  ~7#~=% o&ͦoͦܐԩͣ}!!"8~#(}:$= +*!Z!*B!!:(=2!Z: <2!!!:O::O:!*B! !45(!.+/ 0y0( d!kZ!{Z͈͈o&  :(y ͠|( *"x2y( >28!?"9!!>2 :D]SXN]D [ (!e}̈́A8Q0G: x@!\w# (   yV. V!h6# (*(.(!8}(*(̈́w#>?> w#a{ |͒}͛Ɛ'@'7||}>"C"6# ""͐ͩ*B"[R5*"^#V#^#V#N#FO/o&9O/o&9!9(> (G!9 w#E͊w}8uRB0 >R@RR!+ͨ z R!+ͨ z <!+ͨ z <!+ͨ z <!#ͨ z <!+ͨ z T]KB!z> S>))0 = |JJDMgo>jB0 7?= H\<z5+)+<z {0Gɯgo||H}||/g}/o#}o&K[xAJSJDM!b"!6J"DM'ͬͬdͬ ͬ} wͦWͧ _}8(8J`9{T]=o`9y w >uJ u` }>(; xQ }} ˸T}ٕ(0D=C ,= ( [ 0%D , 7 ͏ ?(8u x O - ; 8˸x X ,-xG}; }م 9; .>#n0[ D = - nx P ,-(-˸G,-; }ٕ? 9.>͏ 8u ?= u+-(>O 0u O 8͏ ?x P , 78ƀ8ƀ8ox٨!دoGOW_gɷɷ|لg{ً_zيWyىOxوG|ٔg{ٛ_zٚWyٙOx٘Gxٸyٹzٺ{ٻ|ټx٨ xx(ͼ ?}ٽÏ }ց; <(; 7D = |٤g{٣_z٢Wy١Ox٠GD u J }x>uu}ƀ/ƀo; -J }0W-J W,}l˸ͨ 8 ; ` x( -ͨ 8J -ͨ 8,J }l8;*!` ! >u` ` u--- J ,,,-xGg?+2n*8t z~,->uxua}.; OJ , ; !U >,k- o&0%,` }g; }؉}颋.:}8c~I$I~L*kٷx˸; }0G,͙<},-(-J ! >0 a` o8 Oþ >um.`1pF,t6|!wS<.z}[|%FXc~ur1}Oٯx(<˸ͨ 8; !~Jͨ 0O!><ͨ 8 =  7 <` O ; 7 0 W-J OT0 j oD,:j !I}袋.}8c~I$I~L!>u` ` 77 ` = O nf^VNF!DLT\I!!53!r1!\!> x #-= o˸xO(- }(x>8(C ,C `iM!>u|; |J>| )=|(DMbo˸ͦ88ͦx(0 8> Mx(>-Ͳ{(ay(Ͱͦ \z(>.Ͳ (Ͱ ~ͦ{>EͲ>+|(|Dg>-Ͳ|/ 0:p# ~# +>0w#,-  60#J˸}րogM| .(C = ~> x0w#xG%P %P ZJDM%P = _~65i+~hìx-Sx9?+{Η@}|C C gZJDM0D ,7}o˸  #yO!@9i&   # w# /w# w#!9! E9!!9~(+F͊!"9!(#>2*Ͳ"|>" :( ͆ *6#w*6#6 !\$![ (̈́( #:~CONTRMKBDLSTAUXUSR>2$*#~ Ͷ$*:> >w###6  #6++p>2S-$Ͷ:*6###ww#w$w#w: ##N#F*B> w#w#[s#r>2S$Ͷ$*6 #-Nw#Fwq#p#6#w#w#w* :( ͒: *^ F* < >26"~͟*-w#ww#͟"~ <@*Ͳ!\  <ʮ!\$> >2*|>! * \$\<(!: [1Á\!(f"> 2:!<"F( #~#6e>!["N>!~8>O6*"w (=(&("( :(N 8y(~#x+% (6*#~[*#~ *~(h#"b=  8 J= B== ͯ}8= ͵}/ͭ !*###~-_~(4Q6*>2>*##w:>*##~*#~(E[ ( ( ( !][ ( ( ((w#(6!]~-#8~>7  [>OkͼMs #rkͼpX á[ [ (( #w(q*#~[ (  *##~6͜O$*#~(08ʦ=ʦ==ʩ=ʬò+###~-_q46͡> *:4^q}Ò*|(M|( M6-#͐ͦ[R8 (G> ͒C~͒#*ͦC!h !lTRUEFALSEͦ!9^#(~#(G~͒#> ͒> Ò "F![(#RR0*4#4> RR *4 #4(>>2$*V(/˖:(#~+ x y2!͵( =( X:(R*:(###~-_-͌X> :("͟"*^˞*V˖0 SRѷR8A* N#F#s#r$ 0})jS\*###w* N#FB ͟r+s> !T]>)j)0 0= UR!#U*^#V#N#F#^#V>">!2DM"~x(L* :O(o:" C}=( ?*-N#Fp+qq#p! * F+N+++V+^Bq#p>>> SRѷR* s#r$ s#r"S"! N#FB(^x * 6#[<(H*! Kq#p##K[! *! 4 #4! x *$ *>w""{_!"*nf}(HR0nf" ^VMDnfutqp*s#r*s#r"* 5KB!>u~#fo{_"*R0RnfR0KqputsrNF( ^VNF^V*SutKqp R*R(~w~wnf ut"6# * *!""*NFy(* "*B0Cnf* [R*"*RS[s#r^#VS>O"w2x2!"" @*>2"!"""!\Ͳ*: !~6go(\R*s#r_2x( s x(T]DMR0 -a%}̈́o*!~6o&͠|ͣ}%^C User break1:% I/O% Run-time% error ͒%, PC=[R"͍% Program aborted*1!͍!Z"ͲI"">"@!")*)n&"%*%!͛E !"%!*@s#r!"'*'*%͇*@^#V!ͯ}oEʞ!!+Ͳ*'*%͇*)*'+n&! N}oEʕ *'!"'X *'*%͇*)*'+n&! a}o!+ m!ͯ}oE!!+ *)*'+n&e.!+Ͳ*'!"'Õ *'*%͇*)*'+n&! a}oEK!*'!"'!!+ m!͛Eʛ!*@^#V!*@s#r*>*@^#V+)))!+ ( "" ""*"!͛*"*͇}oE!!B*"+))) * Ͳ!}2$"!}2$*$&!"*!!͟!EH"![! !q͐b*!" "!!B")) * Ͳ!}2$!!}2$ { Copyright (C) 1984 D.M. Fritz-Rohner } function SS ( var a: vertex ): real ; { } { Name: SS - normalized least squares extremum. } { } { Version: 1.0/Turbo Date: 1984 May 16 } { } { Purpose: Compute normalized root sum square of } { difference between specified and computed } { dependent variable. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { } { SS extremum function, sum of squares } { a[1..m] parameter vector } { ff(a[1..m],x) parameter defined function } { } var i : index ; s : real ; { intermediate variable } begin if n > 0 then begin s := 0.0 ; for i := 1 to n do s := s + sqr(ys[i] - ff(a,xs[i])) ; SS := sqrt (s / n) end else begin writeln(err,'Unexpected Data Order: n < 1') ; quit end end ;  { Copyright (C) 1984 D.M. Fritz-Rohner } function SS ( var a: vertex ): real ; { } { Name: SS - mini-max extremum function. } { } { Version: 1.0/Turbo (tm) Date: 1984 May 16 } { } { Purpose: Compute maximum of absolute value of } { difference between specified and computed } { independent variable values. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { } { SS extremum function, sum of squares } { a[1..m] parameter vector } { ff(a[1..m],x) parameter defined function } { } var i : index ; s : real ; { intermediate variable } sx : real ; begin if n > 0 then begin sx := abs(ys[1] - ff(a,xs[1])) ; if n > 1 then for i := 2 to n do begin s := abs(ys[i] - ff(a,xs[i])) ; if s > sx then sx := s end ; SS := sx end else begin writeln(err,'Unexpected Data Order: n < 1') ; quit end end ;  { Copyright (C) 1984 D.M. Fritz-Rohner } function SS ( var a: vertex ): real ; { } { Name: SS - Pass without change computed function } { } { Version: 1.0/Turbo Date: 1984 May 16 } { } { Purpose: Support application of SIMPLEX procedure } { to functions where no 'fit' criterion is } { defined } { } { This routine simply passes the value defined by the } { function ff. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } begin SS := ff(a,0.0) end ;  { Copyright (C) 1984 D.M. Fritz-Rohner } function ffTest ( var a: vertex ) : boolean ; { } { Name: ffTest - feasibility test } { } { Version: 1.2/TURBO (tm) Date: 1984 July 11 } { } { Purpose: Feasibility test always true. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { ffTest feasibility test result } { a[1..m] parameter vector } { } begin ffTest := true end ;  { Copyright (C) 1984 D.M. Fritz-Rohner } function ffTest ( var a: vertex ) : boolean ; { } { Name: ffTest - feasibility test } { } { Version: 1.2/TURBO (tm) Date: 1984 July 11 } { } { Purpose: Feasibility test for multi-linear fit. } { } { Description: Return true if breakpoint abcissae are } { ordered; } { } { x < x < ... < x } { 1 2 m/2 } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { ffTest feasibility test result } { a[1..m] parameter vector, where; } { } { a[1..m] = [x , y , x , y , ... , x , y ] } { 1 1 2 2 m/2 m/2 } { } var i : integer ; begin if m > 3 then begin i := m - 1 ; while (i > 1) and (a[i-2] < a[i]) do i := i - 2 ; ffTest := i = 1 end else begin writeln(err,'Unreasonable Parameter Vector Order') ; quit end end ; function ff ( var a: vertex ; x : real ) : real ; { } { Name: ff - compute fitting function. } { } { Version: 1.0/TURBO (tm) Date: 1984 May 16 } { } { Purpose: Compute 'Michaelis-Menten' function value. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { Reference: 'Fitting Curves to Data' } { Marco S. Caceci and William P. Cacheris } { BYTE Magazine, Volume 9, Number 4 } { May 1984, pp. 340-362 } { } { } { ff dependent variable, estimated } { a[1..m] parameter vector } { x independent variable } { } begin ff := a[1] * x / (x + a[2]) end ; 200 5e-3 2 0.2 1.0 1.0 3.0 1.0 1.0 6 0.172 1.68 0.250 3.33 0.286 5.00 0.303 6.67 0.334 10.0 0.384 20.0 { Copyright (C) 1984 D.M. Fritz-Rohner } function ff ( var a: vertex ; x : real ) : real ; { } { Name: ff - compute fitting function value. } { } { Version: 1.2/TURBO (tm) Date: 1984 July 11 } { } { Purpose: Compute simple odd polynomial value. } { } { ff = x * (a[1] + x2 * (a[2] .. + x2 * a[m]))..) } { } { where: } { } { x2 = x * x } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { ff dependent variable, estimated } { a[1..m] parameter vector } { x independent variable } { } var x2 : real ; j : integer ; sum : real ; begin x2 := x * x ; { initialize } sum := 0.0 ; for j := 1 to m do sum := a[m-j+1] + x2 * sum ; ff := x * sum end ; 100 1.e-4 3 1.0 1.0 1.0 -.167 1.0 1.0 .00833 1.0 1.0 50 3.141076E-02 2.000000E-02 6.279052E-02 4.000000E-02 9.410831E-02 6.000000E-02 1.253332E-01 8.000000E-02 1.564345E-01 1.000000E-01 1.873813E-01 1.200000E-01 2.181432E-01 1.400000E-01 2.486899E-01 1.600000E-01 2.789911E-01 1.800000E-01 3.090170E-01 2.000000E-01 3.387379E-01 2.200000E-01 3.681246E-01 2.400000E-01 3.971479E-01 2.600000E-01 4.257793E-01 2.800000E-01 4.539905E-01 3.000000E-01 4.817537E-01 3.200000E-01 5.090414E-01 3.400000E-01 5.358268E-01 3.600000E-01 5.620834E-01 3.800000E-01 5.877853E-01 4.000000E-01 6.129071E-01 4.200000E-01 6.374240E-01 4.400000E-01 6.613119E-01 4.600000E-01 6.845471E-01 4.800000E-01 7.071068E-01 5.000000E-01 7.289686E-01 5.200000E-01 7.501111E-01 5.400000E-01 7.705132E-01 5.600000E-01 7.901550E-01 5.800000E-01 8.090170E-01 6.000000E-01 8.270806E-01 6.200000E-01 8.443279E-01 6.400000E-01 8.607420E-01 6.600000E-01 8.763067E-01 6.800000E-01 8.910065E-01 7.000000E-01 9.048271E-01 7.200000E-01 9.177546E-01 7.400000E-01 9.297765E-01 7.600000E-01 9.408808E-01 7.800000E-01 9.510565E-01 8.000000E-01 9.602937E-01 8.200000E-01 9.685832E-01 8.400000E-01 9.759168E-01 8.600000E-01 9.822873E-01 8.800000E-01 9.876883E-01 9.000000E-01 9.921147E-01 9.200000E-01 9.955620E-01 9.400000E-01 9.980267E-01 9.600000E-01 9.995066E-01 9.800000E-01 1.000000E+00 1.000000E+00  { Copyright (C) 1984 D.M. Fritz-Rohner } function ff ( var a: vertex ; x : real ) : real ; { } { Name: ff - compute fitting function. } { } { Version: 1.2/TURBO (tm) Date: 1984 July 11 } { } { Purpose: Compute multi-linear fit. } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { ff dependent variable, estimated } { a[1..m] parameter vector } { x independent variable } { } var i : integer ; x0,y0,x1,y1 : real ; begin if m > 3 then begin i := m - 1 ; while (i > 0) and (x < a[i]) do i := i - 2 ; if i = m - 1 then ff := a[m] else if i < 1 then ff := a[2] else begin x0 := a[i] ; y0 := a[i+1] ; x1 := a[i+2] ; y1 := a[i+3] ; ff := y0 + (y1-y0)*(x-x0)/(x1-x0) end end else begin writeln(err,'Unexpected Parameter Vector Order') ; quit end end ; 200 1e-4 6 0.00000 0.1 1.0 0.00000 0.1 1.0 0.30000 0.1 1.0 0.29552 0.1 1.0 0.60000 0.1 1.0 0.56464 0.1 1.0 7 0.000000 0.0 0.099833 0.1 0.198669 0.2 0.295520 0.3 0.389418 0.4 0.479426 0.5 0.564642 0.6 { Copyright (C) 1984 D.M. Fritz-Rohner } function ff ( var a: vertex ; x : real ) : real ; { } { Name: RF - Root Finder. } { } { Version: 1.1/TURBO (tm) Date: 1984 May 16 } { } { Purpose: Compute absolute value sum for polynomial } { whose root(s) are being sought. Here, the } { roles of some of the program variables are } { exchanged; } { } { ff extremum function } { } { a[1..m] become the roots, } { xs[1..n] become the coefficients, and } { ys[1..n] dummied for consistency } { } { m } { _____ } { \ } { ff = \ abs [ (..(x * a + x ) * a + ..) + x ] } { / 1 i 2 i n } { /____ } { i = 1 } { } { Author: D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } var i,j : integer ; sum : real ; num : real ; begin sum := 0.0 ; for i := 1 to m do begin num := 0.0 ; for j := 1 to n do num := num * a[i] + xs[j] ; { Horner's Rule } sum := sum + abs(num) end ; ff := sum end ; 100 1e-5 2 0.0 1.0 1e-5 0.0 1.0 1e-5 3 0.0 1.0 0.0 -1.0 0.0 -2.0