anchor.diff:          print *,' '
   c
           call ajusp(pa,pp,ax,ay,az,xmat,vect,p,cp0,a,b,c,d
  -     &             ,tetav0,rho0,itp,n1,n2,n3)
  +     &             ,tetav0,rho0,itp,n1,n2,n3,n3p1)
           if(ibreak.eq.1)then
             testmax=testmax*testmax
             go to 2
  
positive.code: c**** 1/RHO*d(RHO*VP*QP)/dz
  c**** AT (I+1/2,J+1/2,K+1/2)
  c
        function derz(q,vp,rho0,n1,n2,n3)
  c
        dimension q(n1,n2,n3),vp(n1,n2,n3),rho0(40)
  c
        common/grid/hx,hy,hz,nx,ny,nz,i2d
        common/coord/i,j,k
  
negative.code: c**** THIS ROUTINE WRITES THE MATRIX AND VECTOR FOR PRESSURE RETRIEVAL
  c
        subroutine ajusp(pa,pp,ax,ay,az,x,v,p,cp0,a,b,c,d
       &                 ,tetav0,rho0,ip,n1,n2,n3,n3p1)
  c
        real ax(n1,n2,n3),ay(n1,n2,n3),az(n1,n2,n3)
       &    ,pa(n1,n2,n3),pp(n1,n2,n3)
       &    ,x(27,n1,n2,n3),v(n1,n2,n3),p(n1,n2,n3)
       &    ,a(n1,n2,n3),b(n1,n2,n3),c(n1,n2,n3),d(n1,n2,n3)
       &    ,tetav0(n3p1),rho0(n3p1)
        integer ip(n1,n2,n3,7)
  c
        common/grid/hx,hy,hz,nx,ny,nz,i2d
  