Sources Esope | Cast3M

compac
C COMPAC    SOURCE    PV        22/04/15    17:10:49     11344          
      SUBROUTINE COMPAC(VAL,NBPAR,KIVPO,KIVLO,NVALL,IPPVV,
     >                  IZROSF,NA,PREC,imasq,iprel,iderl)
******
*  COMPACTE UNE LIGNE DE LA MATRICE. TOUTES SEQUENCES DE IZROFS VALEURS
*  INFERIEURES A PREC SONT ELIMINEES
*
*  Modif janvier 2015 toutes les inconnues de la ligne ont le même ivpo
*  de plus, le premier terme de la ligne est forcement gardé
*
*
******
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
-INC CCHOLE
      DIMENSION VAL(*),KIVPO(*),KIVLO(*),IPPVV(*),imasq(*)
*  d'abord calcul des kivpo kivlo initiaux
*  on suppose que les lignes ont la meme colonne de depart. C'est la responsibilite de assem2
*     write (6,*) ' compac na izrosf ',na,izrosf
*     write (6,*) ' ippvv ',(ippvv(i),i=1,na+1)
      incdeb=iprel-ippvv(2)+1
      incdebi=incdeb
      incpos=incdeb
      incfin=incdeb-1
      ildeb=1
      ilfin=0
      nbpar=0
*  decalage colonne adresse dans la ligne
      idec=-incdebi+1
* on balaye en colonne a partis de incdeb
   1  continue
*  recherche nouveau debut
*  recherche acceleree en utilisant imasq
 
1020  continue
* recherche du saut pour chaque ligne
* verif si premier masque nul
      incdeba=iprel+na
      do 1100 il=1,na
       ipr=(incdeb+ippvv(il)+idec)/masdim+1
       incbas=min(ipr*masdim,iprel+ippvv(il)+idec)
       do inc=incdeb+ippvv(il)+idec,incbas
        if (abs(val(inc)).gt.prec.or.inc.eq.1) goto 1954
       enddo
 
1047  continue
      im=ipr
      do 1045 im=ipr,(iprel-1+ippvv(il)+idec)/masdim+1
**        write (6,*) ' im imasq ',im,imasq(im),imasq(im+1)
          if (imasq(im+1).gt.0) goto 1046
          if (imasq(im+1).lt.0) then
            ipr=-(imasq(im+1)/masdim)+1
**           write (6,*) 'compac acceleration im ipr',im,ipr
            goto 1047
          endif
1045  continue
1046  continue
**    write (6,*) 'compac il incdeb indecba im',il,incdeb,
**   >      im*masdim-idec-ippvv(il),im
      incdeba=min(incdeba,im*masdim-idec-ippvv(il))
1100  continue
**    if (incdeb.ne.incdeba)
**   >  write (6,*)'incdeb incdeba na',incdeb,incdeba,na,ipr,im,ippvv(2)
      incdeb=incdeba
1954  continue
      do 40 incpos=incdeb,iprel-1
       do 45 il=1,na
        inc=incpos-incdebi+ippvv(il)+1
        if (abs(val(inc)).lt.prec.and.inc.ne.1) goto 45
        goto 50
  45   continue
  40  continue
*      write (6,*) ' pas trouve debut '
*  pas trouve de debut
      incdeb=iprel
      goto 200
  50  continue
      incdeb=incpos
      nbpar=nbpar+1
      kivlo(nbpar)=ildeb
      kivpo(nbpar)=incdeb-incdebi+1
* recherche de la fin
* acceleration en utilisant itmasq
      izro=0
      do 20 incpos=incdeb+1,iprel-1
       do 25 il=1,na
        inc=incpos-incdebi+ippvv(il)+1
         if (inc.eq.1.or.abs(val(inc)).gt.prec) then
          izro=0
          goto 26
         endif
  25   continue
       izro=izro+1
  26   continue
       if (izro.ge.izrosf) goto 30
  20  continue
      incpos=iprel-1
  30  continue
      incfin=min(iprel-1,incpos-izro)
**    incfin=min(iprel-1,incpos-1)
      ilfin=ildeb+incfin-incdeb
*     write (6,*) ' incdeb incfin izro iprel',incdeb,incfin,izro,iprel
      ildeb=ilfin+1
      incdeb=incpos+1
      if (incdeb.le.iprel-1) goto 1
 
 200  continue
 
**
*  on a construit kivpo et kivlo il ne reste plus qu'a recopier les valeurs.
*  on a tout fait avant iprel. Il faut donc rajouter la partie triangulaire
**
        nbpar=nbpar+1
        kivlo(nbpar)=ilfin+1
        kivpo(nbpar)=incdeb-incdebi+1
        nvall=kivlo(nbpar)-1
        nvallg=kivpo(nbpar)-1
*       write (6,*) ' compac na nbpar initial ',na,nbpar,nvall,incfin
 
 
*  puis recopier l'ensemble en completant le triangle
 
        ippvv(1)=1
        do 100 il=2,na
         ippvv(il)=(il-1)*nbpar+1
         do 110 nbp=1,nbpar
          kivlo(nbp+(il-1)*nbpar)=kivlo(nbp+(il-2)*nbpar)+nvall+il-1
          kivpo(nbp+(il-1)*nbpar)=kivpo(nbp+(il-2)*nbpar)+nvallg+il-1
 110     continue
 
 100    continue
         nbpar=nbpar*na
         nbpar=nbpar+1
         ippvv(na+1)=nbpar
         kivpo(nbpar)=kivpo(nbpar-1)+na
         kivlo(nbpar)=kivlo(nbpar-1)+na
         nvall=kivlo(nbpar)-1
 
*     write (6,*) ' compac iprel na nbpar nvall ',iprel,na,nbpar,nvall
*     write (6,*) 'nouveau ippvv',(ippvv(i),i=1,na+1)
*     write (6,*) 'kivpo',(kivpo(i),i=1,nbpar)
*     write (6,*) 'kivlo',(kivlo(i),i=1,nbpar)
 
*  deplacement des valeurs
*  fait maintenant dans le programme appelant pour eviter un double mouvement
*     do 300 nbp=1,nbpar-1
*      do 310 iv=kivlo(nbp),kivlo(nbp+1)-1
*         val(iv)=val(iv-kivlo(nbp)+kivpo(nbp))
*310   continue
*300  continue
 
      return
      end