#include "dims.h"
      SUBROUTINE STENMOD(IMX,JMX,C,PHIHM,P,COLATC)
      implicit none
C     ****
C     ****     MODIFY STENCIL TO SET POTENTIAL TO HEELIS VALUE
C     ****       WITHIN AURORAL CIRCLE
C     ****
#include "params.h"
#include "consdyn.h"
!
! Args:
      integer,intent(in) :: imx,jmx
      real,intent(out) :: C(IMX,JMX,*)
      real,intent(in)  :: PHIHM(IMX0,JMX0),COLATC(IMX0,JMX0),
     |                    P(IMX0,JMX0)
!
! Local:
      integer :: nint,i0,j0,i,j,n,jj
C     ****
C     ****     COMPUTE SEPARATION OF GRID POINTS FOR THIS RESOLUTION
C     ****
      NINT = (IMX0-1)/(IMX-1)
      I0 = 1-NINT
      J0 = 1-NINT
C     ****
C     ****     IF (NINT.EQ.1), WE ARE AT HIGHEST RESOLUTION.
C     ****       CORRECT RHS, WHICH IS IN C(10)
C     ****
      IF(NINT.EQ.1)THEN
        DO 3 J = 1,JMX
          DO 3 I = 1,IMX
            C(I,J,10) = P(I,J)*C(I,J,10)+(1.-P(I,J))*C(I,J,9)*
     1        (DLATM/(10.*DTR))**2*PHIHM(I,J)
    3   CONTINUE
      ENDIF
C     ****
C     ****     MODIFY STENCIL, C(I,J,N),N=1,9)
C     ****
      DO 1 J = 1,JMX
        JJ = J0+J*NINT
        DO 2 N = 1,8
          DO 2 I = 1,IMX
            C(I,J,N) = C(I,J,N)*P(I0+I*NINT,JJ)
    2   CONTINUE
        DO 4 I = 1,IMX
          C(I,J,9) = C(I,J,9)*P(I0+I*NINT,JJ)+(1.-P(I0+I*NINT,JJ))*
     1      C(I,J,9)*(DLATM*FLOAT(NINT)/(10.*DTR))**2
    4 CONTINUE
    1 CONTINUE
      RETURN
      END
C
C------------------------------------------------------------------------
      SUBROUTINE STENMD(IMX,JMX,C,PHIHM,P,COLATC)
      
! am_02/02: stenmd: calculate stencil for modified mudpack solver (isolve = 2 in dyncal) 
      
      implicit none
C     ****
C     ****     MODIFY STENCIL TO SET POTENTIAL TO HEELIS VALUE
C     ****       WITHIN AURORAL CIRCLE
C     ****
#include "params.h"
#include "consdyn.h"
      
! globale  unmodified coefficient for using modified mudpack
      real :: cofum
      common/mudmd/cofum(imx0,jmx0,9) 
!
! Args:
      integer,intent(in) :: imx,jmx
      real,intent(out) :: C(IMX,JMX,*)
      real,intent(in)  :: PHIHM(IMX0,JMX0),COLATC(IMX0,JMX0),
     |                    P(IMX0,JMX0)
!
! Local:
      integer :: nint,i0,j0,i,j,n,jj
C     ****
C     ****     COMPUTE SEPARATION OF GRID POINTS FOR THIS RESOLUTION
C     ****
      NINT = (IMX0-1)/(IMX-1)
      I0 = 1-NINT
      J0 = 1-NINT
C     ****
C     ****     IF (NINT.EQ.1), WE ARE AT HIGHEST RESOLUTION.
C     ****       CORRECT RHS, WHICH IS IN C(10)
C     ****
      IF(NINT.EQ.1)THEN
        DO 3 J = 1,JMX
          DO 3 I = 1,IMX
            C(I,J,10) = P(I,J)*C(I,J,10)+(1.-P(I,J))*C(I,J,9)*
     1        (DLATM/(10.*DTR))**2*PHIHM(I,J)
    3   CONTINUE
      ENDIF
C     ****
C     ****     MODIFY STENCIL, C(I,J,N),N=1,9)
C     ****
      if(nint.eq.1) then
        DO 11 J = 1,JMX
          JJ = J0+J*NINT
          DO 21 N = 1,8
            DO 21 I = 1,IMX
              C(I,J,N) = C(I,J,N)*P(I0+I*NINT,JJ)
              cofum(I,J,N) = cofum(I,J,N)*P(I0+I*NINT,JJ)
   21     CONTINUE
          DO 41 I = 1,IMX
            C(I,J,9) = C(I,J,9)*P(I0+I*NINT,JJ)+(1.-P(I0+I*NINT,JJ))*
     1        C(I,J,9)*(DLATM*FLOAT(NINT)/(10.*DTR))**2
            cofum(I,J,9) =cofum(I,J,9)*P(I0+I*NINT,JJ)+
     1                    (1.-P(I0+I*NINT,JJ))*cofum(I,J,9)*
     1                    (DLATM*FLOAT(NINT)/(10.*DTR))**2
   41    CONTINUE
   11    CONTINUE
       else
        DO 1 J = 1,JMX
          JJ = J0+J*NINT
          DO 2 N = 1,8
            DO 2 I = 1,IMX
              C(I,J,N) = C(I,J,N)*P(I0+I*NINT,JJ)
    2     CONTINUE
          DO 4 I = 1,IMX
            C(I,J,9) = C(I,J,9)*P(I0+I*NINT,JJ)+(1.-P(I0+I*NINT,JJ))*
     1        C(I,J,9)*(DLATM*FLOAT(NINT)/(10.*DTR))**2
    4   CONTINUE
    1   CONTINUE
      endif
      
      RETURN
      END
C