2_results (732181), страница 4

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 4)

writeln('***********************************************************');

end;

procedure done; {final message}

begin

Sound(222); Delay(111); Sound(444); Delay(111); NoSound; {beep}

write('* Task processing time '); clock; saveTime;

writeln('* Status: Inverse problem has been successfully evaluated.');

end;

Begin

about;

loadData;

initParameters;

directTask;

for m:=1 to mLast do

begin

if fMulti then

begin

mCur:=m;

clockMessage;

end;

doMinimization; {main part of work}

setApproximationData( Rg, mCur ); {set new approx. data}

resultMessage;

if(( fMulti )AND( m < nPoints ))then reduceSILimits;

end;

done;

End.

П1.5 Модуль глобальных данных EData

Unit EData;

Interface

Uses DOS;

Const

maxPAR = 40; {nodes of approximation max number}

maxFUN = 20; {excitation frequencies max number}

maxSPC = 4; {support approximation values number}

iterImax = 50; {max number of internal iterations}

Const

apSpline = 1; {approximation type identifiers}

apHypTg = 3;

apExp = 2;

apPWCon = 4;

apPWLin = 8;

Type

Parameters = array[ 1..maxPAR ] of real; {Si,Mu data}

Functionals = array[ 1..maxFUN ] of real; {Voltage data}

SpecialPar = array[ 1..maxSPC ] of real; {Special data}

Var

hThick : real; {thickness of slab}

nPoints : integer; {nodes of approximation number within [ 0,H ]}

nLayers : integer; {number of piecewise constant SI layers within[ 0,H ]}

nFreqs : integer; {number of excitation frequencies}

nStab : integer; {required number of true digits in SI estimation}

epsU : real; {relative error of measurement Uvn*}

nApprox : byte; {approximation type identifier}

incVal : real; {step for numerical differ.}

parMaxH : integer; {max number of integration steps}

parMaxX : real; {upper bound of integration}

parEps : real; {error of integration}

derivType: byte; {1 for right; 2 for central}

Var

freqs : Functionals; {frequencies of excitment Uvn*}

Umr, Umi : Functionals; { Re(Uvn*),Im(Uvn*) for "measured" data}

Uer, Uei : Functionals; { Re(Uvn*),Im(Uvn*) for estimated data}

mu : Parameters; {relative permeability nodal values}

si, si0 : Parameters; {conductivity approximation nodal values}

siMin, siMax : Parameters; {conductivity nodal values restrictions}

par0 : SpecialPar; {alfa,si2,beta,gamma - for exp&HypTg}

parMin,parMax: SpecialPar; {-||- min/max}

zLayer : Parameters; {relative borders of slab layers [0,1]}

Var

aG : real; {scale factor for SImin/max}

Ft : real; {current discrepancy functional value}

fMulti : boolean; {TRUE if it isn't an EXP-approximation}

fHypTg : boolean; {TRUE for Hyperbolic tg approximation}

parEqlB : boolean; {TRUE if b[i]=const}

mCur : integer; {current number of approximation nodes}

inFileName : string; {data file name}

outFileName : string; {results file name}

Var

Rg : Parameters; {current SI estimation}

RgS : Parameters; {previous SI estimation}

Gr : array [ 1..2 ,1..maxPAR ] of real; {SI max/min}

Fh : array [ 1..maxPAR , 1..maxFUN ] of real; {current discrepancies}

Type

TTime=record

H, M, S, S100 : word; {hour,min,sec,sec/100}

end;

Var

clk1, clk2 : TTime; {start&finish time}

procedure loadData; {load all user-defined data from file}

Implementation

procedure loadData;

var

i,eqlB : integer;

FF : text;

begin

assign( FF, outFileName ); {clear output file}

rewrite( FF );

close( FF );

assign( FF, inFileName ); {read input file}

reset( FF );

readln( FF );

readln( FF );

readln( FF, hThick, nPoints, nLayers, nFreqs, nStab, epsU, aG, nApprox );

readln( FF );

for i:=1 to nFreqs do read( FF, freqs[i] );

readln( FF );

readln( FF );

readln( FF );

for i:=1 to nPoints do readln( FF, si[i], siMin[i], siMax[i] );

readln( FF );

readln( FF );

readln( FF , incVal, parMaxH, parMaxX, parEps, derivType, eqlB );

readln( FF );

readln( FF );

for i:=1 to maxSPC do readln( FF, par0[i] , parMin[i] , parMax[i] );

readln( FF );

if ( eqlB=0 )then

begin

for i:=1 to nLayers+1 do read( FF, zLayer[i] );

parEqlB:=false;

end

else parEqlB:=true;

close( FF );

for i:=1 to maxPAR do mu[i]:=1;

end;

Var

str : string;

Begin

if( ParamCount = 1 )then str:=ParamStr(1)

else

begin

write('Enter I/O file name, please: ');

readln( str );

end;

inFileName :=str+'.txt';

outFileName:=str+'.lst';

End.

П1.6 Модуль работы с файлами EFile

Unit EFile;

Interface

Uses

DOS, EData;

function isStable( ns : integer; var RG1,RG2 ) : boolean;

function saveResults( ns,iter : integer ) : boolean;

procedure saveExpResults;

procedure saveHypTgResults;

procedure clock;

procedure saveTime;

Implementation

Var

FF : text;

i : byte;

function decimalDegree( n:integer ) : real;{10^n}

var

s:real; i:byte;

begin

s:=1;

for i:=1 to n do s:=s*10;

decimalDegree:=s;

end;

function isStable( ns:integer ; var RG1,RG2 ) : boolean;

var

m : real;

R1 : Parameters absolute RG1;

R2 : Parameters absolute RG2;

begin

isStable:=TRUE;

m:=decimalDegree( ns-1 );

for i:=1 to mCur do

begin

if NOT(( ABS( R2[i]-R1[i] )*m )<=ABS( R2[i]) ) then isStable:=FALSE;

RgS[i]:=Rg[i];

end;

end;

function saveResults( ns , iter : integer ) : boolean;

var

sum : real;

begin

sum:=0;

for i:=1 to nFreqs do sum:=sum + Fh[1,i];

sum:=SQRT( sum/nFreqs );

assign( FF , outFileName );

append( FF );

write( iter:2, ' ', sum:10:7, ' Rg=' );

write( FF , iter:2, ' ', sum:10:7, ' Rg=');

for i:=1 to mCur do

begin

write( Rg[i]:6:3, ' ');

write( FF , Rg[i]:6:3, ' ');

end;

writeln;

writeln( FF );

close( FF );

saveResults:=isStable( ns , Rgs , Rg );

end;

procedure saveExpResults;

begin

assign( FF , outFileName );

append( FF );

writeln( ' siE=',Rg[2]:6:3,' siI=',Rg[1]:6:3,' alfa=',Rg[3]:6:3);

writeln( FF , ' siE=',Rg[2]:6:3,' siI=',Rg[1]:6:3,' alfa=',Rg[3]:6:3);

write( ' SI: ');

write( FF , ' SI: ');

for i:=1 to nPoints do

begin

write( si[i]:6:3,' ');

write( FF , si[i]:6:3,' ');

end;

writeln;

writeln( FF );

close( FF );

end;

procedure saveHypTgResults;

begin

assign( FF , outFileName );

append( FF );

writeln( ' si1=',Rg[2]:6:3,' si2=',Rg[1]:6:3,' beta=',Rg[3]:6:3,' gamma=',Rg[4]:6:3);

writeln( FF , ' si1=',Rg[2]:6:3,' si2=',Rg[1]:6:3,' beta=',Rg[3]:6:3,' gamma=',Rg[4]:6:3);

write( ' SI: ');

write( FF , ' SI: ');

for i:=1 to nPoints do

begin

write( si[i]:6:3,' ');

write( FF , si[i]:6:3,' ');

end;

writeln;

writeln( FF );

close( FF );

end;

procedure clock; {t2 = t2-t1}

var

H1,M1,S1,H2,M2,S2,sec1,sec2 : longint;

begin

GetTime( clk2.H, clk2.M, clk2.S, clk2.S100 ); {current time}

H2:=clk2.H; M2:=clk2.M; S2:=clk2.S; H1:=clk1.H; M1:=clk1.M; S1:=clk1.S;

sec2:= ( H2*60 + M2 )*60 + S2;

sec1:= ( H1*60 + M1 )*60 + S1;

if( sec2 < sec1 )then sec2:=sec2 + 85020; {+23.59.59}

sec2:=sec2 - sec1;

clk2.H := sec2 div 3600; sec2:=sec2 - clk2.H*3600;

clk2.M := sec2 div 60; sec2:=sec2 - clk2.M*60;

clk2.S := sec2;

writeln( clk2.H:2, ':', clk2.M:2, ':', clk2.S:2 );

end;

procedure saveTime;

begin

assign( FF , outFileName );

append( FF );

write( FF ,'* Processing time ',clk2.H:2, ':', clk2.M:2, ':', clk2.S:2 );

close( FF );

end;

End.

П1.7 Модуль решения прямой задачи ВТК для НВТП EDirect

{****************************************************************************}

{ ERIN submodule : EDirect , 15.02.99, (C) 1999 by Nikita U.Dolgov }

{****************************************************************************}

{ Estimates Uvn* for Eddy current testing of inhomogeneous multilayer slab }

{ with surface( flat ) probe. }

{ It can do it using one of five types of conductivity approximation : }

{Spline, Exponential, Piecewise constant, Piecewise linear,Hyperbolic tangent}

{****************************************************************************}

{$F+}

Unit EDirect;

Interface

Uses EData, EMath;

Type

siFunc = function( x:real ) : real;

Var

getSiFunction : siFunc; {for external getting SI estimate}

procedure initConst( par1,par2:integer; par3,par4:real; par5:boolean );

procedure getVoltage( freq : real ; var ur,ui : real ); { Uvn* = ur + j*ui }

procedure setApproximationType( approx : byte ); { type of approx. }

procedure setApproximationItem( SIG:real ; N : byte ); { set SIGMA[ N ]}

procedure setApproximationData( var SIG; nVal : byte ); { SIGMA[1..nVal] }

procedure getApproximationData( var SIG ; var N : byte ); { get SIGMA[ N ]}

Implementation

Const

PI23 = 2000*pi; {2*pi*KHz}

mu0 = 4*pi*1E-7; {magnetic const}

Var

appSigma : Parameters; {conductivity approximation data buffer}

appCount : byte; {size of conductivity approximation data buffer}

appType : byte; {conductivity approximation type identifier}

Type

commonInfo=record

w : real; {cyclical excitation frecuency}

R : real; {equivalent radius of probe}

H : real; {generalized lift-off of probe}

Kr : real; {parameter of probe}

eps : real; {error of integration}

xMax : real; {upper bound of integration}

steps : integer; {current number of integration steps}

maxsteps: integer; {max number of integration steps}

Nlay : integer; {number of layers in slab}

sigma : Parameters; {conductivity of layers}

m : Parameters; {relative permeability of layers}

b : Parameters; {thickness of layers}

zCentre : Parameters; {centre of layer}

end;

procFunc = procedure( x:real; var result:complex);

Var

siB, siC, siD : Parameters; {support for Spline approx.}

cInfo : commonInfo; {one-way access low level info}

function siSpline( x:real ) : real;{Spline approximation}

begin

if( appCount = 1 )then

siSpline := appSigma[ 1 ]

else

siSpline:=Seval( appCount, x, appSigma, siB, siC, siD);

end;

function siExp( x:real ) : real;{Exponential approximation}

begin

siExp:=(appSigma[2]-appSigma[1])*EXP( -appSigma[3]*(1-x) ) + appSigma[1];

end;

function siPWConst( x:real ) : real;{Piecewise constant approximation}

var

dx, dh : real; i : byte;

begin

if( appCount = 1 )then siPWConst := appSigma[ 1 ]

else

begin

dh:=1/( appCount-1 );

dx:=dh/2;

i:=1;

while( x > dx ) do

begin

i:=i + 1;

dx:=dx + dh;

end;

siPWConst:=appSigma[ i ];

end;

end;

function siPWLinear( x:real ) : real;{Piecewise linear approximation}

var

dx, dh : real;

i : byte;

begin

if( appCount = 1 )then siPWLinear := appSigma[ 1 ]

else

begin

dh:=1/( appCount-1 );

dx:=0;

i:=1;

repeat

i:=i + 1;

dx:=dx + dh;

until( x <= dx );

siPWLinear:=appSigma[i-1]+( appSigma[i]-appSigma[i-1] )*( x/dh+2-i);

end;

end;

function siHyperTg( x:real ) : real;{Hyperbolic tangent approximation}

begin

siHyperTg:=appSigma[2]+(appSigma[1]-appSigma[2])*(1+th((appSigma[3]-x)/appSigma[4]))/2;

end;

procedure setApproximationType( approx : byte );

begin

appType := approx;

write('* conductivity approximation type : ');

case approx of

apSpline : begin

writeln('SPLINE');

getSiFunction := siSpline;

end;

apExp : begin

writeln('EXP');

getSiFunction := siExp;

end;

apPWCon : begin

writeln('PIECEWISE CONST');

getSiFunction := siPWConst;

end;

apPWLin : begin

writeln('PIECEWISE LINEAR');

getSiFunction := siPWLinear;

end;

apHypTg : begin

writeln('HYPERBOLIC TANGENT');

getSiFunction := siHyperTg;

end;

end;

end;

procedure setApproximationData( var SIG ; nVal : byte );

var

Sigma : Parameters absolute SIG; i:byte;

begin

appCount := nVal;

for i:=1 to nVal do appSigma[ i ]:=Sigma[ i ];

if( appType = apSpline )then Spline( appCount, appSigma, siB, siC, siD);

end;

procedure getApproximationData( var SIG ; var N : byte );

var

Sigma : Parameters absolute SIG; i:byte;

begin

N := appCount;

for i:=1 to appCount do Sigma[ i ]:=appSigma[ i ];

end;

procedure setApproximationItem( SIG:real ; N : byte );

begin

appSigma[ N ] := SIG;

if( appType = apSpline )then Spline( appCount, appSigma, siB, siC, siD);

end;

procedure functionFi( x:real ; var result:complex );{get boundary conditions function value}

var

beta : array[ 1..maxPAR ]of real;

q : array[ 1..maxPAR ]of complex;

fi : array[ 0..maxPAR ]of complex;

th , z1 , z2 , z3 , z4 , z5 , z6 , z7 : complex;

i : byte;

begin

mkComp( 0, 0, fi[0] );

with cInfo do

for i:=1 to Nlay do

begin

beta[i]:=R*sqrt( w*mu0*sigma[i] ); {calculation of beta}

mkComp( sqr(x), sqr(beta[i])*m[i], z7 ); {calculation of q, z7=q^2}

SqrtC( z7, q[i] );

mulComp( q[i], b[i], z6 ); {calculation of th,z6=q*b}

tanH( z6, th ); {th=tanH(q*b)}

mkComp( sqr(m[i])*sqr(x), 0, z6 ); {z6=m2x2}

SubC( z6, z7, z5); {z5=m2x2-q2}

AddC( z6, z7, z4); {z4=m2x2+q2}

MulC( z5, th, z1); {z1=z5*th}

MulC( z4, th, z2); {z2=z4*th}

mulComp( q[i], 2*x*m[i], z3 ); {z3=2xmq}

SubC( z2, z3, z4 );

MulC( z4, fi[i-1], z5 );

SubC( z1, z5, z6 ); {z6=high}

AddC( z2, z3, z4 );

MulC( z1, fi[i-1], z5 );

SubC( z4, z5, th ); {th=low}

DivC( z6, th, fi[i] );

end;

eqlComp( result, fi[ cInfo.Nlay ] );

end;

procedure funcSimple( x:real; var result:complex );{intergrand function value}

var

z : complex;

begin

with cInfo do

begin

functionFi( x, result );

mulComp( result, exp( -x*H ), z );

mulComp( z, J1( x*Kr ), result );

mulComp( result, J1( x/Kr ), z );

eqlComp( result, z );

end;

end;

procedure funcMax( x:real; var result:complex );{max value; When Fi[Nlay]=1}

var

z1, z2 : complex;

begin

with cInfo do

begin

mkComp( 1,0,z1 );

mulComp( z1, exp(-x*H), z2 );

mulComp( z2, J1( x*Kr ), z1 );

mulComp( z1, J1( x/Kr ), result );

end;

end;

procedure integralBS( func:procFunc ; var result:complex );{integral by Simpson}

var

z , y , tmp : complex;

hh : real;

i, iLast : word;

begin

with cInfo do

begin

hh:=xMax/steps;

iLast:=steps div 2;

nullComp(tmp);

func( 0, z );

eqlComp( y, z );

for i:=1 to iLast do

begin

func( ( 2*i-1 )*hh , z );

deltaComp( tmp, z );

end;

mulComp( tmp, 4, z );

deltaComp( y, z );

nullComp( tmp );

iLast:=iLast-1;

for i:=1 to iLast do

begin

func( 2*i*hh ,z );

deltaComp( tmp, z );

end;

mulComp( tmp, 2, z );

deltaComp( y, z );

func( xmax, z );

deltaComp(y,z);

mulComp( y, hh/3, z );

eqlComp( result, z );

end;

end;{I = h/3 * [ F0 + 4*sum(F2k-1) + 2*sum(F2k) + Fn ]}

procedure integral( F:procFunc; var result:complex );{integral with given error}

var

e , e15 : real;

flag : boolean;

delta , integ1H , integ2H : complex;

begin

with cInfo do

begin

e15 :=eps*15;{ | I2h-I1h |/(2^k -1 )< Eps ; k=4 for Simpson method}

steps:=20;

flag :=false;

integralBS( F, integ2H );

repeat

eqlComp( integ1H, integ2H );

steps:=steps*2;

integralBS( F, integ2H );

SubC( integ2H, integ1H, delta );

e:=Leng( delta );

if( e

until( ( flag ) OR (steps>maxsteps) );

if( flag )then

begin

eqlComp( result, integ2H );

end

else

begin

writeln('Error: Too big number of integration steps.');

halt(1);

end;

end;

end;

procedure initConst( par1, par2 : integer; par3, par4 : real; par5:boolean );

var

i : byte;

bThick, dl, x : real;

const

Ri=0.02; hi=0.005; { radius and lift-off of excitation coil}

Rm=0.02; hm=0.005; { radius and lift-off of measuring coil}

begin

with cInfo do

begin

Nlay :=par1;

xMax :=par3;

maxsteps:=par2;

R :=sqrt( Ri*Rm );

H :=( hi+hm )/R;

Kr :=sqrt( Ri/Rm );

eps :=par4;

bThick :=hThick*0.002/R; {2*b/R [m]}

for i:=1 to Nlay do m[i]:= mu[i];

if par5 then

begin

bThick:=bThick/NLay;

for i:=1 to Nlay do b[i]:=bThick;

dl:=1/NLay;

x:=dl/2; {x grows up from bottom of slab to the top}

for i:=1 to Nlay do

begin

zCentre[i]:=x;

x:=x + dl;

end;

end

else

for i:=1 to Nlay do

begin

b[i]:=( zLayer[i+1]-zLayer[i] )*bThick;

zCentre[i]:=( zLayer[i+1]+zLayer[i] )/2;

end;

end;

end;

procedure init( f:real );{get current approach of conductivity values}

var

i : byte;

begin

with cInfo do

begin

w:=PI23*f;

for i:=1 to Nlay do sigma[i]:=getSiFunction( zCentre[i] )*1E6;

end;

end;

procedure getVoltage( freq : real ; var ur,ui : real );

var

U, U0, Uvn, tmp : complex;

begin

init( freq );

integral( funcSimple, U ); { U =Uvn }

integral( funcMax , U0 ); { U0=Uvn max }

divComp( U, Leng(U0), Uvn ); { Uvn=U/|U0| }

mkComp( 0, 1, tmp ); { tmp=( 0+j1 ) }

MulC( tmp, Uvn, U ); { U= j*Uvn = Uvn* }

ur := U.re;

ui := U.im;

end;

END.

П1.8 Модуль решения обратной задачи ВТК для НВТП EMinimum

Unit EMinimum;

INTERFACE

Uses EData, Crt, EFile, EDirect;

procedure doMinimization;

IMPLEMENTATION

procedure getFunctional( Reg : byte );

var

ur, ui, dur, dui, Rgt : real;

ur2, ui2: Functionals;

i, j, k : byte;

begin

getApproximationData( si , k );

setApproximationData( Rg, mCur );

case Reg of

0 : for i:=1 to nFreqs do {get functional F value}

begin

getVoltage( freqs[i], ur, ui );

Uer[ i ]:=ur; {we need it for dU}

Uei[ i ]:=ui;

Fh[1,i] := SQR( ur-Umr[i] ) + SQR( ui-Umi[i] );

end;

{Right:U'(i)= (U(i+1)-U(i))/h}

1 : for i:=1 to mCur do {get dF/dSI[i] value}

begin

Rgt:=Rg[i]*( 1+incVal ); {si[i]=si[i]+dsi[i]}

setApproximationItem( Rgt, i ); {set new si[i] value}

for j:=1 to nFreqs do

begin {get dUr/dSI,dUi/dSI}

getVoltage( freqs[ j ], ur, ui );

dur:=( ur-Uer[j] )/( Rg[i]*incVal );

dui:=( ui-Uei[j] )/( Rg[i]*incVal );

Fh[i,j]:=2*(dur*(Uer[j]-Umr[j])+dui*(Uei[j]-Umi[j]));

end;

setApproximationItem( Rg[i], i ); {restore si[i] value}

end;

{Central:U'(i)= (U(i+1)-U(i-1))/2h}

2 : for i:=1 to mCur do {get dF/dSI[i] value}

begin

Rgt:=Rg[i]*( 1+incVal ); {si[i]=si[i]+dsi[i]}

setApproximationItem( Rgt, i ); {set new si[i] value}

for j:=1 to nFreqs do getVoltage( freqs[j],ur2[j],ui2[j] );

Rgt:=Rg[i]*( 1-incVal ); {si[i]=si[i]-dsi[i]}

setApproximationItem( Rgt, i ); {set new si[i] value}

for j:=1 to nFreqs do

begin {get dUr/dSI,dUi/dSI}

getVoltage( freqs[ j ], ur, ui );

dur:=( ur2[j]-ur )/( 2*Rg[i]*incVal );

dui:=( ui2[j]-ui )/( 2*Rg[i]*incVal );

Fh[i,j]:=2*(dur*(Uer[j]-Umr[j])+dui*(Uei[j]-Umi[j]));

end;

setApproximationItem( Rg[i], i ); {restore si[i] value}

end;

end;

setApproximationData( si , k );

end;

procedure doMinimization;

const

mp1Max = maxPAR + 1;

mp2Max = maxPAR + 2;

m2Max = 2*( maxPAR + maxFUN );

m21Max = m2Max + 1;

n2Max = 2*maxFUN;

m1Max = maxPAR + n2Max;

n1Max = n2Max + 1;

mm1Max = maxPAR + n1Max;

minDh : real = 0.001; {criterion of an exit from golden section method}

var

A : array [ 1 .. m1Max , 1 .. m21Max ] of real;

B : array [ 1 .. m1Max] of real;

Sx: array [ 1 .. m21Max] of real;

Zt: array [ 1 .. maxPAR] of real;

Nb: array [ 1 .. m1Max] of integer;

N0: array [ 1 .. m21Max] of integer;

a1, a2, dh, r, tt, tp, tl, cv, cv1, cl, cp : real;

n2, n1, mp1, mp2, mm1, m1, m2, m21 : integer;

ain : real;

apn : real;

iq : integer;

k1 : integer;

n11 : integer;

ip : integer;

iterI : integer;

i,j,k : integer;

label

102 ,103 ,104 ,105 ,106 ,107 ,108;

begin

n2:=2*nFreqs; n1:=n2+1; m1:=mCur+n2;

mp1:=mCur+1; mp2:=mCur+2; mm1:=mCur+n1;

m2:=2*( mCur+nFreqs ); m21:=m2+1;

for k:=1 to m1Max do

for i:=1 to m21Max do

A[k,i]:=0;

iterI:=0;

102:

iterI:=iterI+1;

getFunctional( 0 );

for i:=1 to nFreqs do b[i]:= -Fh[1,i];

getFunctional( derivType );

for k:=1 to mCur do

begin

Zt[k]:=Rg[k];

for i:=1 to nFreqs do

begin

A[i,k+1]:=Fh[k,i];

A[i+nFreqs,k+1]:=-A[i,k+1];

end;

for i:=1 to nFreqs do B[i]:=B[i]+Rg[k]*A[i,k+1];

end;

for i:=1 to nFreqs do B[i+nFreqs]:=-B[i];

for i:=n1 to m1 do B[i]:=Gr[1,i-n2]-Gr[2,i-n2];

for i:=1 to m1 do

begin

if i<=n2 then

for k:=2 to mp1 do B[i]:=B[i]-A[i,k]*Gr[2,k-1];

A[i,1]:=-1;

if i>n2 then

begin

A[i,1]:=0;

for k:=2 to mp1 do

if i-n2=k-1 then A[i,k]:=1

else A[i,k]:=0;

end;

for k:=mp2 to m21 do

if k-mp1=i then A[i,k]:=1

else A[i,k]:=0;

end;

k1:=1;

for k:=1 to n2 do

if B[k1]>B[k] then k1:=k;

for k:=1 to mp1 do A[k1,k]:=-A[k1,k];

A[k1,mCur+1+k1]:=0;

B[k1]:=-B[k1];

for i:=1 to n2 do

if i<>k1 then

begin

B[i]:=B[i]+B[k1];

for k:=1 to mm1 do A[i,k]:=A[i,k]+A[k1,k];

end;

for i:=mp2 to m21 do

begin

Sx[i]:=B[i-mp1];

Nb[i-mp1]:=i;

end;

for i:=1 to mp1 do Sx[i]:=0;

Sx[1]:=B[k1];

Sx[mp1+k1]:=0;

Nb[k1]:=1;

103:

for i:=2 to m21 do N0[i]:=0;

104:

for i:=m21 downto 2 do

if N0[i]=0 then n11:=i;

for k:=2 to m21 do

if ((A[k1,n11]N0[k])) then n11:=k;

if A[k1,n11]<=0 then goto 105;

iq:=0;

for i:=1 to m1 do

if i<>k1 then

begin

if A[i,n11]>0 then

begin

iq:=iq+1;

if iq=1 then

begin

Sx[n11]:=B[i]/A[i,n11]; ip:=i;

end

else

begin

if Sx[n11]>B[i]/A[i,n11] then

begin

Sx[n11]:=B[i]/A[i,n11]; ip:=i;

end;

end;

end

else

if iq=0 then

begin

N0[n11]:=n11;

goto 104;

end;

end;

Sx[Nb[ip]]:=0;

Nb[ip]:=n11;

B[ip]:=B[ip]/A[ip,n11];

apn:=A[ip,n11];

for k:=2 to m21 do A[ip,k]:=A[ip,k]/apn;

for i:=1 to m1 do

if i<>ip then

begin

ain:=A[i,n11];

B[i]:=-B[ip]*ain+B[i];

for j:=1 to m21 do A[i,j]:=-ain*A[ip,j]+A[i,j];

end;

for i:=1 to m1 do Sx[Nb[i]]:=B[i];

goto 103;

105:

for k:=1 to mCur do Sx[k+1]:=Sx[k+1]+Gr[2,k];

a1:=0;

a2:=1.;

dh:=a2-a1;

r:=0.618033;

tl:=a1+r*r*dh;

tp:=a1+r*dh;

j:=1;

108:

if j=1 then tt:=tl else tt:=tp;

106:

for i:=1 to mCur do Rg[i]:=Zt[i]+tt*(Sx[i+1]-Zt[i]);

getFunctional( 0 );

cv:=abs(Fh[1,1]);

if nFreqs>1 then

for k:=2 to nFreqs do

begin

cv1:=abs(Fh[1,k]);

if cv

end;

if (j=1) or (j=3) then cl:=cv

else cp:=cv;

if j=1 then

begin

j:=2;

goto 108;

end;

if dh

if cl>cp then

begin

a1:=tl; dh:=a2-a1; tl:=tp; tp:=a1+r*dh ; tt:=tp; cl:=cp; j:=4;

end

else

begin

a2:=tp; dh:=tp-a1; tp:=tl; tl:=a1+r*r*dh; tt:=tl; cp:=cl; j:=3;

end;

goto 106;

107:

if (iterI < iterImax)AND(NOT saveResults( nStab,iterI )) then goto 102;

end;

End.

Приложение 2 - Удельная электрическая проводимость материалов

Приведем сводку справочных данных согласно[7-9].

Приложение 4 - Abstract

The inverse eddy current problem can be described as the task of reconstructing an unknown distribution of electrical conductivity from eddy-current probe voltage measurements recorded as function of excitation frequency. Conductivity variation may be a result of surface processing with substances like hydrogen and carbon or surface heating.

Mathematical reasons and supporting software for inverse conductivity profiling were developed by us. Inverse problem was solved for layered plane and cylindrical conductors.

Because the inverse problem is nonlinear, we propose using an iterative algorithm which can be formalized as the minimization of an error functional related to the difference between the probe voltages theoretically predicted by the direct problem solving and the measured probe voltages.

Numerical results were obtained for some models of conductivity distribution. It was shown that inverse problem can be solved exactly in case of correct measurements. Good estimation of the true conductivity distribution takes place also for measurement noise about 2 percents but in case of 5 percent error results are worse.

56

Характеристики

Список файлов реферата