Télécharger la plaquette Cast3M


fr en

Télécharger moca.eso

Retour à la liste

Numérotation des lignes : 

moca

  1. C MOCA      SOURCE    CHAT      05/01/13    01:46:52     5004
  2.       subroutine moca
  3.       IMPLICIT INTEGER(I-N)
  4.       IMPLICIT REAL*8(A-H,O-Z)
  5.  
  6. -INC PPARAM
  7. -INC CCOPTIO
  8. -INC SMLREEL
  9.       pointeur mlree4.mlreel, mlree5.mlreel
  10.       character*4 mot(1)
  11.       segment itra
  12.         integer ipoi(m),npara,ibon(m)
  13.       endsegment
  14.       segment itrav
  15.         real*8  a(n,n),b(n),c(n),bb(n,n)
  16.         integer is(n)
  17.       endsegment
  18.       data MOT /'POID'/
  19. * lecture de la valeur de depart des paramatres
  20.       call lirobj ('LISTREEL',mlpar,1,iretou)
  21.       mlreel=mlpar
  22.       segact mlreel
  23.       nparr = prog(/1)
  24. * lecture de la valeur theorique
  25.       call lirobj('LISTREEL',mlreel,1,iretou)
  26.       if(ierr.ne.0) return
  27.       segact mlreel
  28.       npmes=prog(/1)
  29. * lecture de la valeur de depart
  30.       call lirobj('LISTREEL',mlree1,1,iretou)
  31.       if(ierr.ne.0) return
  32.       segact mlree1
  33.       if(mlree1.prog(/1).ne.npmes) then
  34.        write(6,*)' erreur 0 npmes mlree1.prog', npmes,mlree1.prog(/1)
  35.        call erreur (19)
  36.        return
  37.       endif
  38.       m=100
  39.       segini itra
  40. *      write(6,*) ' nb de parametres ',nparr
  41.       do 1 i=1,nparr
  42. * lecture des sensibilité par rapport au n parametres
  43.        call lirobj('LISTREEL',mlree2,1,iretou)
  44.        if(ierr.ne.0) then
  45.         write(6,*)  ' erreur 1'
  46.          call erreur (19)
  47.          return
  48.        endif
  49.        segact mlree2
  50.        do iy=1,mlree2.prog(/1)
  51.         if( mlree2.prog(iy).ne.0.d0) go to 10
  52.        enddo
  53.        go to 1
  54.    10  continue
  55.        ibon(i)=1
  56.        if(mlree2.prog(/1).ne.npmes) then
  57.         write(6,*)' erreur 2  npmes mlree2.prog', npmes,mlree2.prog(/1)
  58.          call erreur (19)
  59.          return
  60.        endif
  61.        npara=npara+1
  62.        ipoi(npara)=mlree2
  63.    1  continue
  64.       if( npara.eq.0) then
  65.         write(6,*)  ' erreur 3  npara ',npara
  66.         call erreur (19)
  67.        return
  68.       endif
  69. * lecture de l'option POIDS
  70.       call lirmot(mot,1,itrou,0)
  71.       if( itrou.eq.0) then
  72.         jg=npmes
  73.          isup=1
  74.         segini mlree5
  75.         do io=1,npmes
  76.           mlree5.prog(io)=1.d0
  77.         enddo
  78.       else
  79. * lecture du poids
  80.         call lirobj('LISTREEL',mlree5,1,iretou)
  81.         if(ierr.ne.0)return
  82.         segact mlree5
  83.         isup=0
  84.       endif
  85.       n = npara
  86.       segini itrav
  87.       do 3 i=1,npmes
  88.       do 4 j=1,npara
  89.        mlree3=ipoi(j)
  90.        x1=mlree3.prog(i)*2. *mlree5.prog(i)*mlree5.prog(i)
  91.        do 5 k=1,npara
  92.         mlree4=ipoi(k)
  93.         xk=mlree4.prog(i)
  94.         a(j,k)=a(j,k)+x1*xk
  95.    5   continue
  96.        b(j) = b(j) + x1 * ( mlreel.prog(i)-mlree1.prog(i))
  97.    4  continue
  98.    3  continue
  99. *
  100. **  appel a l'inversion
  101. *
  102. *      write(6,*) ' avant inversion'
  103. *      write(6,*) ( a(1,it),it=1,npara)
  104. *      write(6,*) ( a(2,it),it=1,npara)
  105. *      write(6,*) ( a(3,it),it=1,npara)
  106.  
  107.       diamax=0.
  108.       do io=1,npara
  109.         if( abs(a(io,io)).gt.diamax) diamax=abs(a(io,io))
  110.       enddo
  111.       eps = diamax / 1.e10
  112.       call inver (A,npara,icrit,c,is,eps)
  113.       if(icrit.ne.0) then
  114.        call erreur(21)
  115.        return
  116.       endif
  117. *      write(6,*) ' apres inversion'
  118. *      write(6,*) ( a(1,it),it=1,npara)
  119. *      write(6,*) ( a(2,it),it=1,npara)
  120. *      write(6,*) ( a(3,it),it=1,npara)
  121.  
  122. *      write(6,*) ' b '
  123. *      write(6,*) ( b(it),it=1,npara)
  124.       if(isup.eq.1) then
  125.        segsup mlree5
  126.       else
  127.        segdes mlree5
  128.       endif
  129.       jg=npara
  130.       segini mlree5
  131.       do 50 io=1,npara
  132.       mlree5.prog(io)=0.d0
  133.       do 6 iu=1,npara
  134.       mlree5.prog(io)=mlree5.prog(io)+a(io,iu)*b(iu)
  135.    6  continue
  136.   50  continue
  137. *      write(6,*) 'valeur trouvées dans moca'
  138. *      write(6,*) (mlree5.prog(io),io=1,npara)
  139.       segdes mlreel,mlree1
  140.       do io=1,npara
  141.       mlree3=ipoi(io)
  142.       segdes mlree3
  143.       enddo
  144.       jg=nparr
  145.       segini mlreel
  146.       ia=1
  147.       mlree2=mlpar
  148.       do ib=1,nparr
  149.       if(ibon(ib).eq.1)then
  150.       prog(ib)=mlree5.prog(ia) + mlree2.prog(ib)
  151.       ia=ia+1
  152.        else
  153.        prog(ib)=mlree2.prog(ib)
  154.       endif
  155.       enddo
  156.       segdes mlreel,mlree2
  157.       segsup mlree5
  158.       call ecrobj('LISTREEL',mlreel)
  159.       segsup itrav,itra
  160.       return
  161.       end
  162.  
  163.  
  164.  
  165.  
  166.  
  167.  
  168.  
  169.  
  170.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales