Ìåíþ
Ïîèñê



ðåôåðàòû ñêà÷àòü Ðåøåíèå îáðàòíîé çàäà÷è âèõðåòîêîâîãî êîíòðîëÿ

     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

          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<e15 )then flag:=true;

          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 );                               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]<A[k1,k]) AND (k<>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<cv1 then cv:=cv1;

        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<MinDh then goto 107;

     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].

Ìàòåðèàë

smin ,[ÌÑì/ì]

smax ,[ÌÑì/ì]

Íåìàãíèòíûå ñòàëè

0.4

1.8

Áðîíçû (ÁðÁ, Áð2, Áð9)

6.8

17

Ëàòóíè (ËÑ59, ËÑ62)

13.5

17.8

Ìàãíèåâûå ñïëàâû (ÌË5-ÌË15)

5.8

18.5

Òèòàíîâûå ñïëàâû (ÎÒ4, ÂÒ3-ÂÒ16)

0.48

2.15

Àëþìèíèåâûå ñïëàâû (Â95, Ä16, Ä19)

15.1

26.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.



Ñòðàíèöû: 1, 2, 3, 4, 5




Íîâîñòè
Ìîè íàñòðîéêè


   ðåôåðàòû ñêà÷àòü  Íàâåðõ  ðåôåðàòû ñêà÷àòü  

© 2009 Âñå ïðàâà çàùèùåíû.