c***********************  problem name: control  ***********************
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine a1xy(x,y,u,ux,uy,rl,itag,values)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            real(kind=rknd), dimension(*) :: values
            character(len=80) :: su
            common /val0/k0,ku,kx,ky,kl
cy
        values(k0)=ux
        values(kx)=1.0e0_rknd
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine a2xy(x,y,u,ux,uy,rl,itag,values)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            real(kind=rknd), dimension(*) :: values
            character(len=80) :: su
            common /val0/k0,ku,kx,ky,kl
cy
        values(k0)=uy
        values(ky)=1.0e0_rknd
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine fxy(x,y,u,ux,uy,rl,itag,values)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            real(kind=rknd), dimension(*) :: values
            character(len=80) :: su
            common /atest2/ix,iy,iu(98),gamma,beta,f0,f1,f2,f3,
     +          c0,c1,c2,c3,bdlw,bdup,dbc,ru(87),su(100)
            common /val0/k0,ku,kx,ky,kl
cy
        q=((f3*u+f2)*u+f1)*u+f0
        qu=(3.0e0_rknd*f3*u+2.0e0_rknd*f2)*u+f1
        p=((c3*u+c2)*u+c1)*u+c0
        pu=(3.0e0_rknd*c3*u+2.0e0_rknd*c2)*u+c1
        values(k0)=-q-rl*p
        values(kl)=-p
        values(ku)=-qu-rl*pu
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine gnxy(x,y,u,rl,itag,values)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            real(kind=rknd), dimension(*) :: values
            character(len=80) :: su
            common /val1/k0,ku,kl
cy
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine gdxy(x,y,rl,itag,values)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            real(kind=rknd), dimension(*) :: values
            character(len=80) :: su
            common /val2/k0,kl,klb,kub,kic,kim,kil
            common /atest2/ix,iy,iu(98),gamma,beta,f0,f1,f2,f3,
     +          c0,c1,c2,c3,bdlw,bdup,dbc,ru(87),su(100)
cy
        values(k0)=dbc
        values(klb)=bdlw
        values(kub)=bdup
        values(kil)=(bdlw+bdup)/2.0e0_rknd
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine p1xy(x,y,u,ux,uy,rl,itag,values)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            real(kind=rknd), dimension(*) :: values
            character(len=80) :: su
            common /val0/k0,ku,kx,ky,kl
            common /atest2/ix,iy,iu(98),gamma,beta,f0,f1,f2,f3,
     +          c0,c1,c2,c3,bdlw,bdup,dbc,ru(87),su(100)
cy
        call uexact(x,y,itag,r,rx,ry,rxx,ryy,rxy)
        values(k0)=(u-r)**2+beta*((ux-rx)**2+(uy-ry)**2)+gamma*rl**2
        values(ku)=2.0e0_rknd*(u-r)
        values(kx)=2.0e0_rknd*(ux-rx)*beta
        values(ky)=2.0e0_rknd*(uy-ry)*beta
        values(kl)=2.0e0_rknd*gamma*rl
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine p2xy(x,y,dx,dy,u,ux,uy,rl,itag,jtag,values)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            real(kind=rknd), dimension(*) :: values
            character(len=80) :: su
            common /val0/k0,ku,kx,ky,kl
            common /atest2/ix,iy,iu(98),gamma,beta,f0,f1,f2,f3,
     +          c0,c1,c2,c3,bdlw,bdup,dbc,ru(87),su(100)
cy
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine qxy(x,y,u,ux,uy,rl,itag,values)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            real(kind=rknd), dimension(*) :: values
            common /val3/kf,kf1,kf2,kad
cy
        call uexact(x,y,itag,r,rx,ry,rxx,ryy,rxy)
        values(kf)=r
        values(kf1)=rx
        values(kf2)=ry
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine sxy(rl,s,itag,values)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            real(kind=rknd), dimension(2,*) :: values
            common /val4/j0,js,jl
cy
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine uexact(x,y,itag,u,ux,uy,uxx,uyy,uxy)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            character(len=80) :: su
            common /atest2/ix,iy,iu(98),gamma,beta,f0,f1,f2,f3,
     +          c0,c1,c2,c3,bdlw,bdup,dbc,ru(87),su(100)
cy
        pi=3.141592653589793e0_rknd
        ax=real(ix,rknd)*pi
        ay=real(iy,rknd)*pi
        sx=sin(ax*x)
        sy=sin(ay*y)
        cx=cos(ax*x)
        cy=cos(ay*y)
        u=sx*sy+dbc
        ux=ax*cx*sy
        uy=ay*sx*cy
        uxx=-ax**2*u
        uyy=-ay**2*u
        uxy=ax*ay*cx*cy
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine usrcmd(vx,vy,sf,itnode,ibndry,ip,rp,sp,iu,ru,su)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            integer(kind=iknd), dimension(5,*) :: itnode
            integer(kind=iknd), dimension(7,*) :: ibndry
            integer(kind=iknd), dimension(100) :: ip,iu
            real(kind=rknd), dimension(*) :: vx,vy
            real(kind=rknd), dimension(2,*) :: sf
            real(kind=rknd), dimension(100) :: rp,ru
            character(len=80), dimension(100) :: sp,su
cy
c       enter input mode
c
        call usrset(iu,ru,su)
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine usrfil(file,len)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            integer(kind=iknd), save :: len0
            character(len=80), save, dimension(100) :: file0
            character(len=80), dimension(*) :: file
cy
            data (file0(i),i=  1, 10)/
     +  'n i= 1,n=ix,    a=ix,t=i',
     1  'n i= 2,n=iy,    a=iy,t=i',
     2  'n i= 1,n=gamma, a=g ,t=r',
     3  'n i= 2,n=beta,  a=b ,t=r',
     4  'n i= 3,n=f0    ,a=f0,t=r',
     5  'n i= 4,n=f1    ,a=f1,t=r',
     6  'n i= 5,n=f2    ,a=f2,t=r',
     7  'n i= 6,n=f3    ,a=f3,t=r',
     8  'n i= 7,n=c0    ,a=c0,t=r',
     9  'n i= 8,n=c1    ,a=c1,t=r'/
            data (file0(i),i= 11, 15)/
     +  'n i= 9,n=c2    ,a=c2,t=r',
     1  'n i=10,n=c3    ,a=c3,t=r',
     2  'n i=11,n=bdlw,  a=bl,t=r',
     3  'n i=12,n=bdup,  a=bu,t=r',
     4  'n i=13,n=dbc,   a=bc,t=r'/
c
            data len0/15/
c
        len=len0
        do i=1,len
            file(i)=file0(i)
        enddo
        return
        end
c-----------------------------------------------------------------------
c
c            piecewise lagrange triangle multi grid package
c
c                  edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
        subroutine gdata(vx,vy,sf,itnode,ibndry,ip,rp,sp,iu,ru,su,sxy)
cx
            use mthdef
            implicit real(kind=rknd) (a-h,o-z)
            implicit integer(kind=iknd) (i-n)
            integer(kind=iknd), dimension(5,*) :: itnode
            integer(kind=iknd), dimension(7,*) :: ibndry
            integer(kind=iknd), dimension(100) :: ip,iu
            integer(kind=iknd), save, dimension(9) :: ix,iy
            integer(kind=iknd), save :: ntf,nvf,nbf,ispd,ifirst
            real(kind=rknd), dimension(*) :: vx,vy
            real(kind=rknd), dimension(2,*) :: sf
            real(kind=rknd), dimension(100) :: rp,ru
            character(len=80), dimension(100) :: sp,su
cy
            external sxy
            data ix/0,0,1,2,2,2,1,0,1/
            data iy/1,2,2,2,1,0,0,0,1/
            data ntf,nvf,nbf,ispd,ifirst/8,9,8,1,1/
c
c           common /atest2/ix,iy,iu(98),gamma,beta,f0,f1,f2,f3,
c    +          c0,c1,c2,c3,bdlw,bdup,dbc,ru(87),su(100)
c
c
        if(ip(41)==1) then
            sp(2)='control'
            sp(1)='control'
            sp(3)='control'
            sp(4)='control'
            sp(5)='control_figxxx.ext'
            sp(6)='control_mpixxx.rw'
            sp(7)='control.jnl'
            sp(9)='control_mpixxx.out'
            iu(1)=5
            iu(2)=5
            ru(1)=1.0e-4_rknd
            ru(2)=0.0e0_rknd
            ru(3)=0.0e0_rknd
            ru(4)=0.0e0_rknd
            ru(5)=0.0e0_rknd
            ru(6)=0.0e0_rknd
            ru(7)=1.0e0_rknd
            ru(8)=0.0e0_rknd
            ru(9)=0.0e0_rknd
            ru(10)=0.0e0_rknd
            ru(11)=1.0e0_rknd
            ru(12)=10.0e0_rknd
            ru(13)=0.0e0_rknd
        endif
c
        ip(1)=ntf
        ip(2)=nvf
        ip(3)=nbf
        ip(5)=max(ip(5),ifirst)
        ip(6)=5
        ip(8)=ispd
        ip(7)=0
        rp(3)=0.1e0_rknd
c**
        ip(20)=4
        ip(21)=4
c**
        do i=1,ntf
            itnode(1,i)=9
            itnode(2,i)=i
            itnode(3,i)=i-1
            itnode(4,i)=0
            itnode(5,i)=i
c
            ibndry(1,i)=i
            ibndry(2,i)=i-1
            ibndry(3,i)=0
            ibndry(4,i)=2
            ibndry(5,i)=0
            ibndry(6,i)=0
            ibndry(7,i)=(i+1)/2
            sf(1,i)=0.0e0_rknd
            sf(2,i)=0.0e0_rknd
        enddo
        itnode(3,1)=8
        ibndry(2,1)=8
c
        do i=1,nvf
            vx(i)=real(ix(i),rknd)/2.0e0_rknd
            vy(i)=real(iy(i),rknd)/2.0e0_rknd
        enddo
c
        return
        end