- 01
- 02
- 03
- 04
- 05
- 06
- 07
- 08
- 09
- 10
- 11
- 12
- 13
- 14
- 15
- 16
- 17
- 18
- 19
- 20
- 21
- 22
- 23
- 24
- 25
- 26
- 27
- 28
- 29
- 30
- 31
- 32
- 33
- 34
- 35
- 36
- 37
- 38
- 39
- 40
- 41
- 42
- 43
- 44
- 45
- 46
- 47
- 48
- 49
- 50
- 51
- 52
- 53
- 54
- 55
- 56
- 57
- 58
- 59
- 60
- 61
- 62
- 63
- 64
- 65
- 66
- 67
- 68
- 69
- 70
- 71
- 72
- 73
- 74
- 75
- 76
- 77
- 78
- 79
- 80
1000:
nume:=0;
for i:=1 to mmes do work^[i]{^}:=work^[i]{^}+hde;
iter:=iter+1;
for j:=1 to 2 do
begin
if j=1 then w:=work^[1]{^};
if j=2 then w:=work^[mmes]{^} ;
ab:=2*sqrt(3.*abs(1.-w));
w1:=ab;
for i:=1 to 10 do
begin
w2:=ab * sqrt(1.+w1);
if abs(w2-w1)<1.e-6 then goto 18;
w1:=w2
end;
18:
wkb:=(1.+w2)/w;
if j=1 then rpsi[1]^:=work^[2]{^}*wkb;
if j=2 then rpsi[mmes]^:=work^[m1]{^}*wkb
end;
b:=rpsi[mmes]^;
for i:=1 to m1 do
begin
npsi:=mmes-i;
a:=(12./work^[npsi]{^}-10.0)-1./b;
rpsi[npsi]^:=a;
b:=a;
if a<=1. then goto 30
end;
30:
mcross:=npsi;
rmcros:=a;
b:=rpsi[1]^;
for i:=2 to mcross do
begin
a:=(12./work^[i]{^}-10.0)-1./b;
rpsi[i]^:=a;
b:=a;
if a<0. then nume:=nume+1
end;
if (iter=1) and (nume<num) then writeln('error in initial Eh');
if nume=num then goto 50;
if (not bisec) then writeln('give left energy El');
if (succes) then goto 60;
if nume>num then eh:=enew;
if nume<num then el:=enew;
goto 80;
50:
m3:=mcross+1;
m4:=mcross-1;
a:=1.-0.5/work^[m3]{^};
a1:=a*(1./rpsi[m3]^-rpsi[mcross]^);
b:=1.-0.5/work^[m4]{^};
b1:=b*(rmcros-1./rpsi[m4]^);
de:=(a1-b1)*work^[mcross]{^};
if de>0. then eh:=enew;
if de<0. then el:=enew;
it:=it+1;
goto 70;
60:
it:=0;
70:
if it>=maxit then goto 100;
80:
eold:=enew;
delta:=(el-eh)*0.5;
enew:=eh+delta;
hde:=hsq12*(enew-eold);
succes:=false;
{writeln ('delta=',delta); }
if abs(delta)>1.e-10 then goto 1000;
{writeln('enew=',enew); }
if nume=num then goto 90;
writeln ('degeneration : num-state=nume-state ');
halt;
90:
{ lip:=true; }
100:
succes:=true;
Follow us!