twodint.f File Reference

Go to the source code of this file.

Functions/Subroutines
subroutine	twodint (T, LSP, IART, XA, YA, ZA, NA, IEXP, IER)

Function/Subroutine Documentation

twodint()

subroutine twodint	(	real*8, dimension(lsp,1)	T,
		integer	LSP,
		integer	IART,
		real*8, dimension(1)	XA,
		real*8, dimension(1)	YA,
		real*8, dimension(1)	ZA,
		integer	NA,
		integer, dimension(2)	IEXP,
		integer	IER
	)

      implicit none
      INTEGER iexp(2),iyu,iyo,ixu,ixo,idx,idy,ll,inpy,iexpx1,iexpxn,
     &  iexpy1,iexpyn,lx,ly,inpx,iart,lsp,ier,nx,ny,l,na,one
      REAL*8 t(lsp,1),xa(1),ya(1),za(1)
      REAL*8 z1(4),z2(4)
C      ENTRY ZWEINT (T,LSP,IART,XA,YA,ZA,NA,IEXP,IER)
      ier = 0
      one=1
      nx = t(1,1)
      ny = (t(1,1)-nx)*1000 + 0.1
C
C TESTING INPUT
C--------------
      IF ((nx-2).lt.0) then
         go to 900
      elseif((nx-2).eq.0) then
         go to 30
      else
         go to 10
      endif
   10 DO 20 l = 3,nx
   20 IF ((t(l,1)-t(l-1,1)) .LE. 0) GO TO 900
   30 IF ((ny-2).lt.0) then
         go to 900
      elseif((ny-2).eq.0) then
         go to 60
      else
         go to 40
      endif
   40 DO 50 l = 3,ny
   50 IF ((t(1,l)-t(1,l-1)) .LE. 0) GO TO 900
   60 IF (na .LE. 0) GO TO 900
C
C DEFINING THE CONTROL VALUES
C---------------------------
      inpx = iart/10
      inpy = iart - inpx*10 + 0.1
      iexpx1 = iexp(1)/10
      iexpxn = iexp(1) - iexpx1*10
      iexpy1 = iexp(2)/10
      iexpyn = iexp(2) - iexpy1*10
      IF (nx-2 .LT. inpx) inpx = nx - 2
      IF (ny-2 .LT. inpy) inpy = ny - 2
      IF (iexpx1 .GT. inpx) iexpx1 = inpx
      IF (iexpxn .GT. inpx) iexpxn = inpx
      IF (iexpy1 .GT. inpy) iexpy1 = inpy
      IF (iexpyn .GT. inpy) iexpyn = inpy
C
C SUCCESSIVE PROCESSING THE INTERPOLATION EXIGENCES
C-------------------------------------------------------
      DO 400 l = 1,na
      lx = 2
C
C SETTING REFERENCE POINTS (LX,LY) 
C---------------------------------
  200 IF (xa(l) .LT. t(lx,1)) GO TO 220
      lx = lx + 1
      IF ((lx-nx).le.0) then
         go to 200
      else
         go to 210
      endif
  210 lx = nx
  220 DO 230 ly = 2,ny
  230 IF (ya(l) .LT. t(1,ly)) GO TO 235
      ly = ny
  235 iyu = ly - inpy
      iyo = ly + inpy - 1
      IF (iyu .GE. 2) GO TO 240
      iyu = 2
      iyo = iyu + inpy
  240 IF (iyo .GT. ny) iyo = ny
      ixu = lx - inpx
      ixo = lx + inpx - 1
      IF (ixu .GE. 2) GO TO 245
      ixu = 2
      ixo = ixu + inpx
  245 IF (ixo .GT. nx) ixo = nx
      idx = ixo - ixu + 1
      IF (ixu .LT. ixo) GO TO 270
      IF (iyu .LT. iyo) GO TO 250
C
C CONSTANT INTERPOLATION
C------------------------
      IF (lx .GT. 2 .AND. xa(l) .LT. t(nx,1)) lx = lx - 1
      IF (ly .GT. 2 .AND. ya(l) .LT. t(1,ny)) ly = ly - 1
      za(l) = t(lx,ly)
      GO TO 400
C
C LINEAR AND  QUADRATIC INTERPOLATION USING ONEDINT (ONEDIMENSIONAL)
C---------------------------------------------------------------------
C
C INTERPOLATION ONLY IN Y-DIRECTION
C
  250 idy = 0
      DO 260 ll = iyu,iyo
      idy = idy + 1
      z1(idy) = t(1,ll)
  260 z2(idy) = t(lx,ll)
      GO TO 300
C
C INTERPOLATION ONLY IN X-DIRECTION
C
  270 IF (iyu .LT. iyo) GO TO 280
      CALL onedint(t(ixu,1),t(ixu,ly),idx,xa(l),za(l),one,inpx,iexp(1),
     1           ier)
      IF (ier.eq.0) then
         go to 400
      else
         go to 900
      endif
C
C 1.INTERPOLATION STEP IN X-DIRECTION
C
  280 idy = 0
      DO 290 ll = iyu,iyo
      idy = idy + 1
      z1(idy) = t(1,ll)
      CALL onedint (t(ixu,1),t(ixu,ll),idx,xa(l),z2(idy),one,inpx,
     1    iexp(1),ier)
      IF (ier.eq.0) then
         go to 290
      else
         go to 900
      endif
  290 CONTINUE
C
C 1.OR 2.INTERPOLATION STEP IN Y-DIRECTION
C
  300 CALL onedint (z1,z2,idy,ya(l),za(l),one,inpy,iexp(2),ier)
      IF (ier.eq.0) then
         go to 400
      else
         go to 900
      endif
C
C RETURN BY NORMAL PROCEEDING
C--------------------------------
  400 CONTINUE
      ier = 0
      RETURN
C
C ERROR RETURN
C-------------
  900 ier = -1
      RETURN