c*********************** problem name: domains ************************
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
common /atest2/iu(100),a1x,a1y,a1u,a2x,a2y,a2u,
+ bux,buy,cu0,cu1,ru(90),su(100)
cy
values(k0)=a1x*ux+a1y*uy+a1u*u
values(ku)=a1u
values(kx)=a1x
values(ky)=a1y
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
common /atest2/iu(100),a1x,a1y,a1u,a2x,a2y,a2u,
+ bux,buy,cu0,cu1,ru(90),su(100)
cy
values(k0)=a2x*ux+a2y*uy+a2u*u
values(ku)=a2u
values(kx)=a2x
values(ky)=a2y
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 /val0/k0,ku,kx,ky,kl
common /atest2/iu(100),a1x,a1y,a1u,a2x,a2y,a2u,
+ bux,buy,cu0,cu1,ru(90),su(100)
cy
values(k0)= - bux*ux - buy*uy - cu0 - cu1*u
values(ku)= - cu1
values(kx) = - bux
values(ky) = - buy
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
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
common /val2/k0,kl,klb,kub,kic,kim,kil
cy
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
common /val0/k0,ku,kx,ky,kl
cy
values(k0)=1.0e0_rknd
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
common /val0/k0,ku,kx,ky,kl
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
s=x**2+y**2
values(kf)=real(itag,rknd)
values(kad)=s
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 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
iu1=iu(1)
c
call usrset(iu,ru,su)
c
ip(41)=0
if(iu1/=iu(1)) ip(41)=-1
c
ispd=ip(8)
a1x=ru(1)
a1y=ru(2)
a1u=ru(3)
a2x=ru(4)
a2y=ru(5)
a2u=ru(6)
bux=ru(7)
buy=ru(8)
cu0=ru(9)
cu1=ru(10)
if(ispd==1) then
if((a1y/=a2x).or.(a1u/=0.0e0_rknd).or.
+ (a2u/=0.0e0_rknd).or.(bux/=0.0e0_rknd).or.
1 (buy/=0.0e0_rknd)) ip(8)=0
else
if((a1y==a2x).and.(a1u==0.0e0_rknd).and.
+ (a2u==0.0e0_rknd).and.(bux==0.0e0_rknd).and.
1 (buy==0.0e0_rknd)) ip(8)=1
endif
if(ispd/=ip(8)) ip(5)=max(ip(5),1)
c
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=domain,a= d,t=i',
1 'n i= 1,n=a1x, a=x1,t=r',
2 'n i= 2,n=a1y, a=y1,t=r',
3 'n i= 3,n=a1u, a=u1,t=r',
4 'n i= 4,n=a2x, a=x2,t=r',
5 'n i= 5,n=a2y, a=y2,t=r',
6 'n i= 6,n=a2u, a=u2,t=r',
7 'n i= 7,n=bux, a=bx,t=r',
8 'n i= 8,n=buy, a=by,t=r',
9 'n i= 9,n=cu0, a=c0,t=r'/
data (file0(i),i= 11, 20)/
+ 'n i=10,n=cu1, a=c1,t=r',
1 's n=domain,v=1,l="texas"',
2 's n=domain,v=2,l="doughnut"',
3 's n=domain,v=3,l="cmos device"',
4 's n=domain,v=4,l="lake superior"',
5 's n=domain,v=5,l="dd20 logo"',
6 's n=domain,v=6,l="at&t logo"',
7 's n=domain,v=7,l="north sea"',
8 's n=domain,v=8,l="airfoil"',
9 's n=domain,v=9,l="planter"'/
data (file0(i),i= 21, 30)/
+ 's n=domain,v=10,l="fan"',
1 's n=domain,v=11,l="arc"',
2 's n=domain,v=12,l="spiral"',
3 's n=domain,v=13,l="ucsd logo"',
4 's n=domain,v=14,l="nozzle"',
5 's n=domain,v=15,l="pltmg"',
6 's n=domain,v=16,l="crack"',
7 's n=domain,v=17,l="monterey bay"',
9 's n=domain,v=18,l="siam logo"',
+ 's n=domain,v=19,l="mexico"'/
data (file0(i),i= 31, 34)/
+ 's n=domain,v=20,l="ellipse"',
1 's n=domain,v=21,l="hole"',
2 's n=domain,v=22,l="british columbia"',
3 's n=domain,v=23,l="israel"'/
c
data len0/34/
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 box(nvf,vx,vy,nbf,ibndry,sf)
cx
use mthdef
implicit real(kind=rknd) (a-h,o-z)
implicit integer(kind=iknd) (i-n)
integer(kind=iknd), dimension(7,*) :: ibndry
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
cy
c rescale domain to lie inside the unit square
c
xmin=vx(1)
xmax=xmin
ymin=vy(1)
ymax=ymin
do i=1,nvf
xmin=min(xmin,vx(i))
xmax=max(xmax,vx(i))
ymin=min(ymin,vy(i))
ymax=max(ymax,vy(i))
enddo
c
c compute scaled coordinates
c
xmid=(xmin+xmax)/2.0e0_rknd
ymid=(ymin+ymax)/2.0e0_rknd
scale=1.0e0_rknd/max(xmax-xmin,ymax-ymin)
do i=1,nvf
vx(i)=0.5e0_rknd+scale*(vx(i)-xmid)
vy(i)=0.5e0_rknd+scale*(vy(i)-ymid)
enddo
do i=1,nbf
if(ibndry(3,i)==0) cycle
sf(1,i)=0.5e0_rknd+scale*(sf(1,i)-xmid)
sf(2,i)=0.5e0_rknd+scale*(sf(2,i)-ymid)
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 :: iprob,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 iprob,ifirst/1,1/
c
if(ip(41)==1) then
sp(4)='domains'
sp(5)='domains_figxxx.ext'
sp(6)='domains_mpixxx.rw'
sp(7)='domains.jnl'
sp(9)='domains_mpixxx.out'
c
iu(1)=15
c
c initialize as laplacian with dirichlet b.c.
c
do i=1,10
ru(i)=0.0e0_rknd
enddo
ru(1)=1.0e0_rknd
ru(5)=1.0e0_rknd
ru(9)=1.0e0_rknd
ip(8)=1
endif
c
ip(5)=max(ip(5),ifirst)
ip(6)=iprob
ip(20)=5
c
if(iu(1)==1) call gd1(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==2) call gd2(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==3) call gd3(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==4) call gd4(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==5) call gd5(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==6) call gd6(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==7) call gd7(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==8) call gd8(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==9) call gd9(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==10) call gd10(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==11) call gd11(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==12) call gd12(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==13) call gd13(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==14) call gd14(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==15) call gd15(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==16) call gd16(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==17) call gd17(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==18) call gd18(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==19) call gd19(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==20) call gd20(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==21) call gd21(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==22) call gd22(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
if(iu(1)==23) call gd23(vx,vy,sf,itnode,
+ ibndry,ip,rp,sp,iu,ru,su,sxy)
c
nvf=ip(2)
nbf=ip(3)
call box(nvf,vx,vy,nbf,ibndry,sf)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd1(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), dimension(12) :: ix,iy
integer(kind=iknd), dimension(5) :: ixc,iyc
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: x,y,xc,yc,hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c domain in the shape of texas.
c dirichlet boundary conditions on all boundary edges
c
external sxy
data ix/2,50,50,90,90,158,159,160,130,90,62,40/
data iy/97,97,160,160,125,121,102,62,0,32,60,60/
data ixc/51,169,106,133,82/
data iyc/54,82,12,42,50/
data ntf,nvf,nbf/1,12,12/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
sp(2)='texas'
sp(1)='texas'
sp(3)='texas'
c
itnode(1,1)=1
itnode(2,1)=12
itnode(3,1)=0
itnode(4,1)=0
itnode(5,1)=1
do i=1,nvf
vx(i)=real(ix(i),rknd)
vy(i)=real(iy(i),rknd)
enddo
do i=1,nbf
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
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)=0
enddo
ibndry(2,nbf)=1
ibndry(3,7)=2
ibndry(3,8)=4
ibndry(3,9)=3
ibndry(3,10)=5
ibndry(3,11)=1
do i=1,nbf
if(ibndry(3,i)>0) then
i1=ibndry(1,i)
i2=ibndry(2,i)
i3=ibndry(3,i)
xc=real(ixc(i3),rknd)
yc=real(iyc(i3),rknd)
call centre(vx(i1),vy(i1),vx(i2),vy(i2),
+ xc,yc,sf(1,i),sf(2,i))
endif
enddo
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c refine long curved edges
c
call sklutl(1_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd2(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 :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c annular region between circles of radius one and six
c dirichlet boundary conditions on inner circle
c neumann boundary conditions on outer circle
c there are two irregular octagons concentric with the
c circles. there are 24 regions, and three groups of
c eight regions are similar.
c
external sxy
data ntf,nvf,nbf/24,32,56/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
sp(2)='doughnut'
sp(1)='doughnut'
sp(3)='doughnut'
c
pi=3.141592653589793e0_rknd
pi4=pi/4.0e0_rknd
c
do i=1,8
r2=2.0e0_rknd
if((i/2)*2==i) r2=4.0e0_rknd
r3=3.0e0_rknd
if((i/2)*2==i) r3=5.0e0_rknd
r4=6.0e0_rknd
ang=real(i-1,rknd)*pi4
k=(i-1)*4+1
vx(k)=cos(ang)
vy(k)=sin(ang)
vx(k+1)=r2*vx(k)
vy(k+1)=r2*vy(k)
vx(k+2)=r3*vx(k)
vy(k+2)=r3*vy(k)
vx(k+3)=r4*vx(k)
vy(k+3)=r4*vy(k)
enddo
do i=1,nbf
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
c
do i=1,32
ibndry(1,i)=i
ibndry(2,i)=i+4
if(i+4>32) ibndry(2,i)=i+4-32
ibndry(3,i)=0
ibndry(4,i)=0
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=0
enddo
do i=1,8
k=(i-1)*4+1
ibndry(3,k)=1
ibndry(4,k)=2
ibndry(3,k+3)=1
ibndry(4,k+3)=1
do j=1,3
k=32+(i-1)*3+j
l=(i-1)*4+j
ibndry(1,k)=l
ibndry(2,k)=l+1
ibndry(3,k)=0
ibndry(4,k)=0
ibndry(5,k)=0
ibndry(6,k)=0
ibndry(7,k)=0
enddo
enddo
c
c comput itnode, find symmetry
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
call sklutl(2_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
c
c label regions
c
itnode(5,1)=3
itnode(5,5)=2
itnode(5,9)=1
do i=1,ntf
j=abs(itnode(3,i))
if(j>0) itnode(5,i)=itnode(5,j)
enddo
c
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd3(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(35) :: b1,b2,b4
integer(kind=iknd), save, dimension(30) :: ix,iy
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c an irregular region composed of six subregions, two quite small,
c which define a typical cmos device
c the data was provided by wolfgang fichtner and donald j. rose
c
external sxy
data ntf,nvf,nbf/6,30,35/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
data (ix(i),i= 1, 30)/
+ 0, 180, 240, 255, 315, 405, 735, 780,
1 825, 1185, 1260, 1275, 1320, 1500, 1500, 1260,
2 1185, 825, 780, 735, 315, 255, 735, 1000,
3 0, 0, 1000, 1500, 0, 1500/
data (iy(i),i= 1, 30)/
+ -11, -12, -31, 25, 37, 41, 40, 21,
1 19, 14, 1, -55, -39, -35, -63, -63,
2 -81, -83, -80, -58, -56, -39, -150, -150,
3 -39, -750, -750, -750,-1000,-1000/
c
data b1/ 1, 2, 3, 4, 5, 6, 7, 8, 9,10,11,12,13,14,15,16,16,
+ 17,18,19,20,21,22,22,25,20,23,24,25,26,27,15,28,29,26/
data b2/ 2, 3, 4, 5, 6, 7, 8, 9,10,11,12,13,14,15,16,12,17,
+ 18,19,20,21,22, 3,25, 1,23,24,27,26,27,28,28,30,30,29/
data b4/ 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 1, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 1, 1, 2, 1/
c
sp(2)='cmos device'
sp(1)='cmos device'
sp(3)='cmos device'
c
do i=1,nvf
vx(i) =real(ix(i),rknd)*2.9563e-4
vy(i) =real(iy(i),rknd)*2.9563e-4
enddo
do i=1,nbf
ibndry(1,i)=b1(i)
ibndry(2,i)=b2(i)
ibndry(3,i)=0
ibndry(4,i)=b4(i)
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=0
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd4(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(303) :: ix,iy
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c region in the shape of lake superior
c there are five subregions and six islands
c the data was provided by r. bruce simpson
c
external sxy
data ntf,nvf,nbf/5,303,313/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
data (ix(i),i= 1, 50)/
+ -891541481,-888470650,-872705746,-867803288,-827638721,
1 -775894356,-733998108,-720648766,-698716116,-691908741,
2 -665943623,-647050762,-639824581,-617760324,-596836376,
3 -605006313,-627510357,-620834494,-635938740,-643797827,
4 -643896866,-633085680,-592166805,-588535309,-553069878,
5 -545711184,-527769947,-468826580,-443255520,-398677874,
6 -362327838,-351568484,-311511421,-260318851,-242016149,
7 -205767846,-194694924,-165644574,-143858552,-118571806,
8 -89657599, -46563160, -7156760, 32215631, 60890228,
9 60928792, 39424559, 28667879, 7170380, 8968640/
data (ix(i),i= 51, 100)/
+ 10765760, -14367931, -32358369, -46769229, -93919420,
1 -104951930,-108823633,-101584518, -94241017, -61523801,
2 -50619459, -28907421, -50666702, -61600220, -39797199,
3 -18084030, 28943431, 54311031, 90588319, 94328558,
4 105342555, 116419089, 131092036, 131172514, 142146897,
5 193172145, 222260547, 233405209, 240735531, 244494963,
6 259090638, 262498355, 273434877, 284545517, 293664694,
7 292019439, 299319220, 331354308, 338635898, 399603295,
8 407055664, 464906311, 490399790, 544955301, 598504639,
9 616538048, 634761524, 663559580, 653090858, 646236992/
data (ix(i),i= 101, 150)/
+ 648030853, 640949774, 641000891, 648826694, 670686483,
1 685469055, 692969561, 736259174, 741950417, 748102570,
2 766219997, 773394203, 769277287, 761990309, 757881975,
3 746499538, 746040964, 766680956, 777094603, 772988129,
4 767073679, 763567352, 762461424, 755321169, 752808619,
5 753453922, 771314287, 786032724, 792287445, 773488998,
6 737138224, 726983166, 718857956, 711526585, 704416895,
7 693369579, 696761417, 706480742, 708971071, 729173613,
8 732301521, 720742178, 709879112, 705625916, 683496714,
9 654260111, 635381174, 664933491, 653851652, 670529938/
data (ix(i),i= 151, 200)/
+ 670097637, 574411726, 556869459, 478860712, 446368170,
1 420757627, 360757804, 338736439, 328043842, 317628479,
2 303819609, 282868814, 285980701, 268455315, 254680014,
3 244376707, 233904195, 219831228, 216162419, 209189844,
4 202351665, 195406818, 157050800, 150021160, 125474477,
5 73096502, 45220131, 17386571, -6945330, -17374951,
6 -20846450, -34761500, -45182380, -55581152, -52159631,
7 -52211881, -48796061, -20926540, -20975170, -52576512,
8 -56044650, -77162850, -94730818, -98336190,-115933979,
9 -122879457,-115743721, -69940150, -59331119, -69731730/
data (ix(i),i= 201, 250)/
+ -97591662,-111755574,-111940694,-129516613,-154830980,
1 -193663633,-186377645,-175828385,-172029650,-154273140,
2 -154171550,-157623363,-199851763,-235065794,-256451440,
3 -253518987,-275629210,-379806066,-433735085,-498349476,
4 -592565298,-703822851,-745170450,-820339012,-872353649,
5 -608656883,-608186007,-579297161,-568086433,-559666157,
6 -566355705,-583187723,-583458996,-579201603,-587263441,
7 -591520166,-602223587,-576770115,-572977018,-549145699,
8 -544981623,-553187752,-553359890,-556836748,-573244572,
9 -534449291,-531934547,-531668663,-529240274,-519116497/
data (ix(i),i= 251, 300)/
+ -521965599,-522571039,-527064848,-259048891,-237567210,
1 -209000087,-184055698,-173409450,-162766552,-155639541,
2 -130752015, -95321161, -88289100, -91850460,-120228684,
3 -116768241,-152397847,-170117927,-191505325,-209371305,
4 -205889106,-188140857,-216642499,-245249653,-252317262,
5 -262850904, -17432991, 3484280, 13937110, 48804179,
6 62737727, 73254800, 82016379, 66322219, 55831802,
7 38390779, 22677940, 3490070, -3493550, -1747930,
8 -10487590, -17459050, 438131142, 437850761, 455353165,
9 483649635, 498189831, 508861971, 516141987, 512747717/
data (ix(i),i= 301, 303)/
+ 498668671, 488142157, 459645033/
data (iy(i),i= 1, 50)/
+ 166159201, 155386531, 143963063, 133102524, 136420846,
1 155212080, 158664119, 174001110, 178389883, 167510176,
2 182378805, 202851081, 202577615, 212360430, 190416360,
3 164209175, 143823195, 127681434, 112335944, 96729130,
4 94082898, 91031909, 116069555, 115944648, 88274997,
5 90689111, 82183212, 106961095, 111578405, 152849865,
6 157335603, 151815784, 161627734, 192534435, 208145809,
7 218255091, 239298439, 260162497, 281135678, 296818018,
8 328396988, 349393415, 365208673, 359945369, 349439454,
9 338853741, 338789296, 341413832, 333451533, 322866678/
data (iy(i),i= 51, 100)/
+ 317574883, 301698995, 285845375, 275293469, 206681108,
1 174987423, 135313761, 132622755, 148459792, 174775326,
2 190618157, 201153159, 174739587, 153603923, 180003631,
3 185259604, 179981732, 166811454, 153733039, 132581246,
4 111473274, 85079402, 69305450, 58719772, 53513449,
5 54008722, 59652221, 43917981, 41371560, 33483499,
6 33695781, 49629170, 49798739, 39389670, 39542049,
7 28925619, 29050219, 71991080, 72132570, 113154662,
8 105386281, 117447603, 115508866, 114505672, 142784142,
9 146074021, 144090629, 150483525, 139480269, 128620267/
data (iy(i),i= 101, 150)/
+ 80988193, 75423539, 70125639, 59818262, 60653287,
1 55932897, 50931108, 63324392, 58265841, 42624259,
2 46068209, 49038869, 59460950, 59138030, 69562429,
3 79670852, 90255553, 117676401, 128750885, 139171314,
4 149511147, 146702886, 130745316, 127777922, 143574929,
5 170116031, 176216900, 171585441, 193089068, 208134031,
6 211820865, 195477509, 216339254, 218681216, 215734005,
7 220580411, 226020551, 252929902, 279541850, 312215924,
8 322956347, 341012549, 343199420, 358924270, 373907781,
9 388636851, 414405775, 479170609, 489328432, 516503143/
data (iy(i),i= 151, 200)/
+ 527087212, 523458147, 517555094, 515094805, 535365105,
1 566463518, 718588781, 765763617, 776130056, 773271894,
2 765061188, 764676332, 785909224, 790898466, 780084085,
3 769334936, 769178009, 776918602, 790102768, 790009642,
4 779334068, 776600647, 773530149, 778755760, 794423151,
5 815268230, 825750160, 830988121, 852151489, 841573811,
6 844224453, 836309528, 838982296, 846954727, 831064415,
7 815185833, 794003391, 783356619, 746306610, 704035950,
8 714633512, 693549252, 688350058, 672492552, 667316341,
9 677953959, 693780470, 741152763, 772867489, 788788366/
data (iy(i),i= 201, 250)/
+ 794227743, 762562370, 736098242, 725644398, 641176081,
1 631007004, 652095509, 651977348, 678403473, 699396944,
2 709982538, 715309381, 699898672, 689787054, 668937254,
3 631832743, 573957253, 501914406, 476767540, 457404757,
4 412823153, 316436315, 275800443, 212452579, 186391926,
5 143142223, 156373394, 173914695, 184134316, 175908267,
6 170830834, 166099274, 158160555, 155366385, 152992165,
7 155788350, 140265095, 195023239, 200193572, 210003328,
8 204571271, 196889234, 191596723, 195948660, 192254853,
9 217481208, 228527403, 236995435, 245395660, 246141505/
data (iy(i),i= 251, 300)/
+ 235634065, 227706409, 222547889, 504862738, 517772865,
1 533265162, 546199608, 548729992, 551267290, 556491137,
2 572152567, 587786102, 582453489, 577180147, 561485386,
3 550874090, 524700975, 524875212, 514524031, 504154873,
4 498817253, 498604679, 491014862, 478181410, 480935097,
5 489040852, 788645267, 799221802, 799227524, 791356945,
6 794052219, 780864763, 772967958, 770248508, 775501299,
7 772803307, 778066111, 772757530, 756878948, 746292973,
8 746296310, 764827394, 445106268, 455691290, 466757298,
9 472859049, 462697220, 463022566, 457950878, 452547693/
data (iy(i),i= 301, 303)/
+ 446819878, 441208935, 440393257/
c
sp(2)='lake superior'
sp(1)='lake superior'
sp(3)='lake superior'
c
do i=1,303
ibndry(1,i)=i-1
ibndry(2,i)=i
ibndry(3,i)=0
ibndry(4,i)=2
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=0
vx(i)=real(ix(i),rknd)*1.0e-7_rknd
vy(i)=real(iy(i),rknd)*1.0e-7_rknd
enddo
ibndry(1,1)=225
ibndry(1,226)=237
ibndry(1,238)=245
ibndry(1,246)=253
ibndry(1,254)=276
ibndry(1,277)=292
ibndry(1,293)=303
do i=304,313
ibndry(3,i)=0
ibndry(4,i)=0
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=0
enddo
ibndry(1,304)=153
ibndry(2,304)=299
ibndry(1,305)=154
ibndry(2,305)=296
ibndry(1,306)=189
ibndry(2,306)=291
ibndry(1,307)=188
ibndry(2,307)=277
ibndry(1,308)=254
ibndry(2,308)=217
ibndry(1,309)=216
ibndry(2,309)=256
ibndry(1,310)=221
ibndry(2,310)=249
ibndry(1,311)=246
ibndry(2,311)=241
ibndry(1,312)=229
ibndry(2,312)=245
ibndry(1,313)=237
ibndry(2,313)=18
c
do i=1,nbf
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd5(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), dimension(200) :: ix,iy,ib1,ib2,ib3,ib4,
+ ib6,isf1,isf2
integer(kind=iknd), dimension(5) :: ixc,iyc
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: x,y,xc,yc,hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c domain in the shape of the dd20 logo
c dirichlet boundary conditions on all boundary edges
c
external sxy
data (ix(i),i=1,100)/
+ 385, 3192, 3516, 5999,10000,10000, 5999, 3516, 3192, 385,
1 1108, 2083, 4301, 5276, 5276, 4301, 2083, 1108, 0, 0,
2 6384, 6384, 6598, 6736, 6736, 6874, 6912, 6930, 6948, 6967,
3 6983, 6996, 7003, 7008, 7010, 7008, 7003, 6996, 6983, 6967,
4 6948, 6930, 6912, 6874, 6754, 6754, 6885, 6939, 6985, 7021,
5 7048, 7065, 7087, 7103, 7116, 7127, 7136, 7143, 7148, 7143,
6 7136, 7127, 7116, 7103, 7087, 7065, 7048, 7021, 6985, 6939,
7 6885, 6598, 7332, 7470, 7470, 7608, 7647, 7665, 7683, 7701,
8 7718, 7730, 7738, 7743, 7745, 7743, 7738, 7730, 7718, 7701,
9 7683, 7665, 7647, 7608, 7488, 7488, 7619, 7674, 7719, 7756/
data (ix(i),i=101,146)/
+ 7783, 7799, 7821, 7838, 7850, 7861, 7870, 7878, 7883, 7878,
1 7870, 7861, 7850, 7838, 7821, 7799, 7783, 7756, 7719, 7674,
2 7619, 7332, 8354, 8476, 8641, 8430, 8630, 8470, 8354, 8238,
3 8078, 8278, 8067, 8232, 9112, 9234, 9399, 9189, 9389, 9229,
4 9112, 8996, 8836, 9036, 8825, 8990/
data (iy(i),i=1,100)/
+ 3379, 1758, 3379, 3379, 3379, 6620, 6620, 6620, 8241, 6620,
1 2516, 1953, 1953, 2516, 7483, 8046, 8046, 7483, 5562, 4437,
2 4437, 5562, 4937, 4937, 5488, 5488, 5477, 5470, 5457, 5437,
3 5413, 5386, 5359, 5322, 5268, 5213, 5177, 5150, 5122, 5099,
4 5079, 5066, 5059, 5048, 5048, 4937, 4937, 4944, 4959, 4975,
5 4995, 5013, 5041, 5068, 5095, 5122, 5159, 5195, 5268, 5341,
6 5377, 5413, 5441, 5468, 5495, 5522, 5541, 5561, 5577, 5591,
7 5599, 5599, 4937, 4937, 5488, 5488, 5477, 5470, 5457, 5437,
8 5413, 5386, 5359, 5322, 5268, 5213, 5177, 5150, 5122, 5099,
9 5079, 5066, 5059, 5048, 5048, 4937, 4937, 4944, 4959, 4975/
data (iy(i),i=101,146)/
+ 4995, 5013, 5041, 5068, 5095, 5122, 5159, 5195, 5268, 5341,
1 5377, 5413, 5441, 5468, 5495, 5522, 5541, 5561, 5577, 5591,
2 5599, 5599, 5170, 4937, 4937, 5275, 5601, 5601, 5382, 5601,
3 5601, 5275, 4937, 4937, 5170, 4937, 4937, 5275, 5601, 5601,
4 5382, 5601, 5601, 5275, 4937, 4937/
data (isf1(i),i=1,24)/
+ 3192, 3192, 0, 0, 0, 0, 0, 3192, 3192, 3192,
1 3192, 0, 3192, 3192, 3192, 3192, 3192, 3192, 3192, 3192,
2 3192, 3192, 3192, 3192/
data (isf2(i),i=1,24)/
+ 5000, 5000, 0, 0, 0, 0, 0, 5000, 5000, 5000,
1 5000, 0, 5000, 5000, 5000, 5000, 5000, 5000, 5000, 5000,
2 5000, 5000, 5000, 5000/
data (ib1(i),i=1,100)/
+ 1, 2, 3, 4, 5, 6, 7, 7, 9, 10,
1 4, 8, 11, 12, 13, 14, 15, 16, 17, 18,
2 19, 20, 21, 22, 23, 24, 25, 26, 27, 28,
3 29, 30, 31, 32, 33, 34, 35, 36, 37, 38,
4 39, 40, 41, 42, 43, 44, 45, 46, 47, 48,
5 49, 50, 51, 52, 53, 54, 55, 56, 57, 58,
6 59, 60, 61, 62, 63, 64, 65, 66, 67, 68,
7 69, 70, 71, 72, 21, 73, 74, 75, 76, 77,
8 78, 79, 80, 81, 82, 83, 84, 85, 86, 87,
9 88, 89, 90, 91, 92, 93, 94, 95, 96, 97/
data (ib1(i),i=101,153)/
+ 98, 99, 100, 101, 102, 103, 104, 105, 106, 107,
1 108, 109, 110, 111, 112, 113, 114, 115, 116, 117,
2 118, 119, 120, 121, 122, 54, 123, 124, 125, 126,
3 127, 128, 129, 130, 131, 132, 133, 134, 104, 135,
4 136, 137, 138, 139, 140, 141, 142, 143, 144, 145,
5 146, 125, 5/
data (ib2(i),i=1,100)/
+ 11, 13, 4, 5, 6, 7, 8, 15, 17, 19,
1 21, 3, 12, 2, 14, 4, 16, 9, 18, 10,
2 20, 1, 22, 7, 24, 25, 26, 27, 28, 29,
3 30, 31, 32, 33, 34, 35, 36, 37, 38, 39,
4 40, 41, 42, 43, 44, 45, 46, 47, 48, 49,
5 50, 51, 52, 53, 54, 55, 56, 57, 58, 59,
6 60, 61, 62, 63, 64, 65, 66, 67, 68, 69,
7 70, 71, 72, 23, 23, 74, 75, 76, 77, 78,
8 79, 80, 81, 82, 83, 84, 85, 86, 87, 88,
9 89, 90, 91, 92, 93, 94, 95, 96, 97, 98/
data (ib2(i),i=101,153)/
+ 99, 100, 101, 102, 103, 104, 105, 106, 107, 108,
1 109, 110, 111, 112, 113, 114, 115, 116, 117, 118,
2 119, 120, 121, 122, 73, 73, 124, 125, 126, 127,
3 128, 129, 130, 131, 132, 133, 134, 123, 133, 136,
4 137, 138, 139, 140, 141, 142, 143, 144, 145, 146,
5 135, 145, 137/
data (ib3(i),i=1,24)/
+ 1,1,0,0,0,0,0,1,1,1,1,0,1,1,1,1,1,1,1,1,1,1,1,1/
data (ib4(i),i=1,153)/
+ 2,2,0,2,2,2,0,2,2,2,0,0,2,2,2,2,2,2,2,2,2,2,0,0,2,2,2,2,2,2,
1 2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,
2 2,2,2,2,2,2,2,2,2,2,2,2,2,2,0,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,
3 2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,
4 2,2,2,2,2,0,2,2,2,2,2,2,2,2,2,2,2,2,0,2,2,2,2,2,2,2,2,2,2,2,
5 2,0,0/
data (ib6(i),i=1,153)/
+ 1, 2, 3, 4, 5, 6, 7, 8, 9,10,11,12, 1, 1, 2, 2, 8, 8, 9, 9,
1 10,10,11,11, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
2 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
3 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1,
4 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
5 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
6 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1,
7 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0/
data ntf,nvf,nbf/3,146,153/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
sp(2)='dd20'
sp(1)='dd20'
sp(3)='dd20'
c
do i=1,nvf
vx(i)=real(ix(i),rknd)
vy(i)=real(iy(i),rknd)
enddo
do i=1,nbf
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
ibndry(1,i)=ib1(i)
ibndry(2,i)=ib2(i)
ibndry(3,i)=0
ibndry(4,i)=ib4(i)
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=ib6(i)
enddo
do i=1,24
sf(1,i)=real(isf1(i),rknd)
sf(2,i)=real(isf2(i),rknd)
ibndry(3,i)=ib3(i)
enddo
c
c uncomment to delete the ddxx symbols
c
cc nbf=24
cc nvf=22
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c refine long curved edges
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
itnode(5,1)=3
itnode(5,2)=2
itnode(5,3)=1
itnode(5,4)=1
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd6(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(159) :: b1,b2,b3,b4,b5
integer(kind=iknd), save, dimension(25) :: c1,c2,c4
integer(kind=iknd), save, dimension(135) :: ix,iy
integer(kind=iknd), dimension(93) :: ixp,iyp
integer(kind=iknd), save :: ntf,ncf,nbf,nvf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save, dimension(93) :: xp,yp
real(kind=rknd), save :: hmax,grade,x,y
character(len=80), dimension(100) :: sp,su
cy
c an irregular region composed of 25 regions, in the shape of,
c the at&t logo.
c the data was provided by r. kent smith
c
c
external sxy
data ntf,nvf,nbf,ncf/25,135,159,93/
c
data (ix(i),i= 1, 50)/
+ -29088581, -2115240, 30000000, 0, -40919900,
1 -31212461, -24240000, -11470000, -5369100, 20187540,
2 27160001, 43160391, 33551569, 21120000, 10151761,
3 -5109150, -13640000, -21440001, -29284019, -36108899,
4 -49742250, -41352282, -32990000, -4100000, 1129340,
5 24291310, 34749999, 50882888, 45282569, 28659999,
6 22388289, 1792850, -8660000, -31520000, -38492460,
7 -46628761, -55230289, -44090900, -37990000, -510000,
8 5761710, 24727149, 35179999, 56014829, 52477179,
9 34890001, 24437151, 5216990, -2100000, -36550000/
data (ix(i),i= 51, 100)/
+ -45265570, -53460388, -58688860, -48302279, -39940000,
1 -1400000, 6444020, 27805979, 35650001, 58897629,
2 56731110, 36570001, 27162440, 5342280, -3020000,
3 -39630001, -46602459, -57378292, -59963450, -45582852,
4 -35130000, -4560000, 3284020, 59990859, 59261298,
5 5425570, -3290000, -37790000, -46505570, -59552770,
6 -59416080, -41115570, -32400000, -7110000, -137540,
7 21172869, 30760000, 59416080, 59990859, 35500000,
8 26784430, -1595980, -9440000, -32820001, -39795649,
9 -59963450, -57063389, -30675910, -11144090, 13666790/
data (ix(i),i= 101, 135)/
+ 21340001, 57063389, 58897629, 28210001, 22629480,
1 -10966530, -34503360, -58897629, -52477179, -26954920,
2 -16485080, 7010050, 10500000, 52477179, 56014829,
3 19870000, 14635080, -15375080, -25844920, -56014829,
4 -45282569, 45282569, 50882888, -50882888, -33551569,
5 33551569, 43160391, -43160391, -30000000, 30000000,
6 0, 0, 400000000, 0,-400000000/
data (iy(i),i= 1, 50)/
+ 52477179, 52805071, 51961522, 60000000, 43881221,
1 44595580, 44899998, 44899998, 44633632, 42354422,
2 42049999, 41679502, 49742250, 49850001, 49179149,
3 47371769, 46849999, 46849999, 47192478, 47918129,
4 33551569, 34211750, 34649999, 34649999, 34421680,
5 32256639, 31800001, 31795161, 39363539, 39600000,
6 39271309, 37097809, 36550000, 36550000, 36854420,
7 37759221, 23443871, 24183631, 24449999, 24449999,
8 24121311, 22247810, 21700001, 21502080, 29088581,
9 29500000, 28952191, 26783471, 26400001, 26400001/
data (iy(i),i= 51, 100)/
+ 26780529, 27239430, 12474700, 13811750, 14250000,
1 14250000, 13907520, 11892480, 11550000, 11448539,
2 19534090, 19349999, 18856970, 16688250, 16250000,
3 16250000, 16554420, 17542300, 2093970, 3702190,
4 4250000, 4250000, 3907520, 1047140, 9386070,
5 6630530, 6250000, 6250000, 6630530, 7312160,
6 -8350390, -6480530, -6100000, -6100000, -6404420,
7 -8331420, -8750000, -8350390, -1047140, -950000,
8 -1330530, -3657520, -4000000, -4000000, -3756410,
9 -2093970, -18541020, -16861030, -16861030, -18732049/
data (iy(i),i= 101, 135)/
+ -19000000, -18541020, -11448539, -11200000, -11394880,
1 -13924609, -13719200, -11448539, -29088581, -28201380,
2 -28201380, -29089079, -29150000, -29088581, -21502080,
3 -21350000, -21441381, -22548621, -22548621, -21502080,
4 -39363539, -39363539, -31795161, -31795161, -49742250,
5 -49742250, -41679502, -41679502, -51961522, -51961522,
6 -60000000, 400000000, 0,-400000000, 0/
c
data (ixp(i),i= 1, 50)/
+ -15328120, 15529140, -15022800, -27729549, -8416640,
1 23670449, 38567259, 15625629, -9366600, -25365739,
2 -38567259, -37176881, -1482840, 29515669, 48231411,
3 25519841, -3426400, -35009551, -48231411, -41043358,
4 2630160, 29946401, 54378471, 29656401, 1563520,
5 -40911942, -54378471, -44126878, 2525740, 31724260,
6 57955551, 31859760, 1166880, -43119550, -58088861,
7 -40363598, -634260, 59771681, 1071940, -42151942,
8 -59815040, -36761940, -3620450, 25961871, 59815040,
9 31138060, -5514260, -36309950, -59771681, -20910001/
data (ixp(i),i= 51, 93)/
+ 17501060, 58088861, 25418041, -22738979, -58088861,
1 -21719999, 8754760, 54378471, 17252140, -20610001,
2 -54378471, 48231411, -48231411, 38567259, -38567259,
3 15529140, -15529140, 31795161, -32678339, 44236641,
4 -43881221, 51697750, -51697750, 56381559, -56381559,
5 59088469, -59177141, 60000000, -60000000, 59177141,
6 -59177141, 56558490, -56558490, 51697750, -51697750,
7 44236641, -44236641, 31795161, -31795161, 282842712,
8 282842712,-282842712,-282842712/
data (iyp(i),i= 1, 50)/
+ 53786950, 57955551, 58088861, 44823861, 44833379,
1 42126141, 45962672, 49682131, 46980562, 46935658,
2 45962672, 34540360, 34592891, 31914210, 35689371,
3 39517770, 36687050, 36626141, 35689371, 24383380,
4 24367771, 21837051, 25357101, 29362950, 26495931,
5 26495180, 25357101, 14140360, 14164340, 11635660,
6 15529140, 19226660, 16359640, 16326140, 15022800,
7 4112950, 4164340, 5229340, 6345180, 6345180,
8 4707550, -6195180, -6176140, -8645300, -4707550,
9 -1045180, -3914340, -3939080, -5229340, -16520000/
data (iyp(i),i= 51, 93)/
+ -18932990, -15022800, -11248730, -14284290, -15022800,
1 -28110001, -29134769, -25357101, -21372850, -22639999,
2 -25357101, -35689371, -35689371, -45962672, -45962672,
3 -57955551, -57955551, 50882888, 50320230, 40535412,
4 40919900, 30452299, 30452299, 20521209, 20521209,
5 10418890, 9902860, 0, 0, -9902860,
6 -9902860, -20028410, -20028410, -30452299, -30452299,
7 -40535412, -40535412, -50882888, -50882888, 282842712,
8 -282842712,-282842712, 282842712/
c
data (b1(i),i=1,100)/
+ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10,
1 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
2 21, 22, 23, 24, 25, 26, 27, 28, 29, 30,
3 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
4 41, 42, 43, 44, 45, 46, 47, 48, 49, 50,
5 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
6 61, 62, 63, 64, 65, 66, 67, 68, 69, 70,
7 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
8 81, 82, 83, 84, 85, 86, 87, 88, 89, 90,
9 91, 92, 93, 94, 95, 96, 97, 98, 99,100/
data (b1(i),i=101,159)/
+ 101,102,103,104,105,106,107,108,109,110,
1 111,112,113,114,115,116,117,118,119,120,
2 121,122,123,124,125,126,127,128,129,130,
3 131, 3, 1, 12, 5, 28, 21, 44, 37, 60,
4 53, 74, 69, 88, 81,102, 97,114,109,122,
5 121,126,125,132,133,134,135,132,134/
data (b2(i),i=1,100)/
+ 2, 3, 4, 1, 6, 7, 8, 9, 10, 11,
1 12, 13, 14, 15, 16, 17, 18, 19, 20, 5,
2 22, 23, 24, 25, 26, 27, 28, 29, 30, 31,
3 32, 33, 34, 35, 36, 21, 38, 39, 40, 41,
4 42, 43, 44, 45, 46, 47, 48, 49, 50, 51,
5 52, 37, 54, 55, 56, 57, 58, 59, 60, 61,
6 62, 63, 64, 65, 66, 67, 68, 53, 70, 71,
7 72, 73, 74, 75, 76, 77, 78, 79, 80, 69,
8 82, 83, 84, 85, 86, 87, 88, 89, 90, 91,
9 92, 93, 94, 95, 96, 81, 98, 99,100,101/
data (b2(i),i=101,159)/
+ 102,103,104,105,106,107,108, 97,110,111,
1 112,113,114,115,116,117,118,119,120,109,
2 122,123,124,121,126,127,128,125,130,131,
3 129, 13, 20, 29, 36, 45, 52, 61, 68, 75,
4 80, 89, 96,103,108,115,120,123,124,127,
5 128,130,129,133,134,135,132, 4,131/
data (b3(i),i=1,100)/
+ 1, 0, 2, 3, 0, 4, 0, 5, 0, 6,
1 0, 7, 0, 8, 0, 9, 0, 10, 0, 11,
2 0, 12, 0, 13, 0, 14, 0, 15, 0, 16,
3 0, 17, 0, 18, 0, 19, 0, 20, 0, 21,
4 0, 22, 0, 23, 0, 24, 0, 25, 0, 26,
5 0, 27, 0, 28, 0, 29, 0, 30, 0, 31,
6 0, 32, 0, 33, 0, 34, 0, 35, 0, 36,
7 0, 37, 0, 38, 0, 39, 0, 40, 0, 41,
8 0, 42, 0, 43, 0, 44, 0, 45, 0, 46,
9 0, 47, 0, 48, 0, 49, 0, 50, 0, 51/
data (b3(i),i=101,159)/
+ 0, 52, 0, 53, 0, 54, 0, 55, 0, 56,
1 0, 57, 0, 58, 0, 59, 0, 60, 0, 61,
2 0, 62, 0, 63, 0, 64, 0, 65, 0, 66,
3 67, 68, 69, 70, 71, 72, 73, 74, 75, 76,
4 77, 78, 79, 80, 81, 82, 83, 84, 85, 86,
5 87, 88, 89, 90, 91, 92, 93, 0, 0/
data (b4(i),i=1,100)/
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
2 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
3 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
4 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
5 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
6 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
7 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
8 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
9 0, 0, 0, 0, 0, 0, 0, 0, 0, 0/
data (b4(i),i=101,159)/
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
2 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
3 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
4 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
5 0, 0, 0, 2, 2, 2, 2, 0, 0/
data (b5(i),i=1,100)/
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
2 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
3 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
4 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
5 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
6 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
7 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
8 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
9 0, 0, 0, 0, 0, 0, 0, 0, 0, 0/
data (b5(i),i=101,159)/
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
2 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
3 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
4 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
5 0, 0, 0, 1, 1, 1, 1, 0, 0/
c
data (c1(i),i=1,25)/
+ 4, 1, 20, 5, 36, 21, 52, 37, 68, 53,
1 80, 69, 96, 81,108, 97,120,109,124,121,
2 128,125,131,134,134/
data (c2(i),i=1,25)/
+ 4,133, 20,135, 36,137, 52,139, 68,141,
1 80,143, 96,145,108,147,120,149,124,151,
2 128,153,130,155,159/
data (c4(i),i=1,25)/
+ 1, 2, 1, 2, 1, 2, 1, 2, 1, 2,
1 1, 2, 1, 2, 1, 2, 1, 2, 1, 2,
2 1, 2, 1, 3, 3/
c
data hmax,grade/0.05e0_rknd,1.5e0_rknd/
c
sp(2)='at&t logo'
sp(1)='at&t logo'
sp(3)='at&t logo'
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
do i=1,nvf
vx(i)=real(ix(i),rknd)*1.0e-7_rknd
vy(i)=real(iy(i),rknd)*1.0e-7_rknd
enddo
do i=1,ncf
xp(i)=real(ixp(i),rknd)*1.0e-7_rknd
yp(i)=real(iyp(i),rknd)*1.0e-7_rknd
enddo
c
do i=1,nbf
ibndry(1,i)=b1(i)
ibndry(2,i)=b2(i)
ibndry(3,i)=b3(i)
ibndry(4,i)=b4(i)
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=b5(i)
c
if(ibndry(3,i)>0) then
i1=ibndry(1,i)
i2=ibndry(2,i)
i3=ibndry(3,i)
call centre(vx(i1),vy(i1),vx(i2),vy(i2),
+ xp(i3),yp(i3),sf(1,i),sf(2,i))
else
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
endif
enddo
c
do i=1,ntf
itnode(1,i)=c1(i)
itnode(2,i)=c2(i)
itnode(3,i)=0
itnode(5,i)=c4(i)
enddo
c
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd7(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(499) :: ix,iy
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c domain in the shape of the north sea.
c the domain contains 29 islands.
c the skeleton has four regions.
c
c data for this domain was provided by gabriel wittum and
c wolfgang hoffman from the university of stuttgart.
c
external sxy
data ntf,nvf,nbf/4,499,531/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
data (ix(i),i= 1, 80)/
+ 8820, 8030, 0, 0, 1260, 3510, 3530, 3810,
1 3840, 3780, 3640, 3700, 3860, 3950, 4140, 4270,
2 4440, 4600, 4680, 4710, 4960, 5040, 5210, 5320,
3 5420, 5530, 5890, 6000, 6030, 6190, 6250, 6520,
4 6550, 6680, 6820, 6770, 6770, 6820, 6880, 6990,
5 7040, 7120, 7180, 7230, 7260, 7210, 7120, 7040,
6 7810, 7810, 7700, 7560, 7480, 7480, 7420, 7530,
7 7560, 7640, 7750, 7840, 7890, 7970, 8050, 8190,
8 8300, 8440, 8520, 8490, 8410, 8600, 8630, 8630,
9 8580, 8520, 8410, 8330, 8300, 8300, 8300, 8360/
data (ix(i),i= 81, 160)/
+ 8410, 8470, 8520, 8520, 8520, 8630, 8710, 8770,
1 8880, 8960, 9070, 9150, 9260, 9320, 9480, 9590,
2 9730, 9840, 9890, 10000, 10000, 9890, 9840, 9780,
3 9730, 9700, 9700, 9700, 9640, 9640, 9590, 9530,
4 9510, 9480, 9530, 9670, 9730, 9780, 9840, 9860,
5 9920, 9950, 9920, 9860, 9860, 9840, 9810, 9780,
6 9780, 9700, 9640, 9560, 9510, 9480, 9450, 9420,
7 9420, 9400, 9400, 9420, 9480, 9480, 9480, 9480,
8 9420, 9370, 9320, 9260, 9230, 9150, 9100, 9010,
9 8990, 8930, 8880, 8790, 8710, 8680, 8660, 8660/
data (ix(i),i= 161, 240)/
+ 8770, 8820, 8930, 8960, 8930, 8880, 8820, 8790,
1 8880, 8960, 9070, 9210, 9290, 9370, 9450, 9530,
2 9620, 9670, 9780, 9840, 9890, 9920, 10000, 10000,
3 9860, 9810, 9730, 9700, 9640, 9640, 9560, 9450,
4 9510, 9480, 9420, 9400, 9320, 9210, 9180, 8990,
5 8990, 9040, 9010, 9010, 8960, 8880, 8880, 8740,
6 8710, 8740, 8740, 8820, 8850, 8850, 8900, 8880,
7 8930, 8930, 8770, 8790, 8820, 9010, 9010, 8930,
8 8850, 8740, 8660, 8660, 8770, 8820, 9010, 8960,
9 8960, 9010, 8900, 8190, 8080, 8050, 8140, 8410/
data (ix(i),i= 241, 320)/
+ 8410, 8300, 8330, 7890, 7970, 7950, 7950, 7840,
1 7890, 7260, 7290, 7210, 7180, 7100, 7100, 7210,
2 6360, 6550, 6330, 6220, 6160, 6250, 5920, 6050,
3 5840, 5750, 5730, 5640, 5620, 5590, 5590, 5670,
4 4900, 4960, 5370, 5420, 5370, 5290, 5180, 5180,
5 5040, 5040, 4900, 4770, 4820, 4680, 4600, 4600,
6 4470, 4440, 4300, 4030, 3970, 3890, 3810, 3950,
7 3150, 3370, 3640, 3510, 3320, 3150, 2960, 2930,
8 3040, 3100, 3040, 2960, 2930, 2900, 2300, 2490,
9 2600, 2660, 2710, 2710, 2490, 2470, 2410, 2490/
data (ix(i),i= 321, 400)/
+ 2330, 2220, 2220, 2030, 2000, 1950, 1920, 1920,
1 1970, 2030, 1530, 1560, 1700, 1750, 1560, 1530,
2 1450, 1450, 1480, 1370, 1370, 1340, 1340, 1260,
3 1230, 6380, 6280, 6270, 6270, 6250, 6330, 6330,
4 9700, 9750, 9860, 9840, 9810, 9750, 9670, 9480,
5 9420, 9400, 9420, 9400, 9420, 9510, 9560, 9640,
6 8580, 8520, 8520, 8580, 8250, 8250, 8300, 8330,
7 8330, 8300, 8250, 8140, 8190, 8270, 8300, 8300,
8 8250, 8190, 8220, 8550, 8580, 8580, 8490, 8470,
9 8380, 8220, 8140, 8190, 8880, 8880, 8820, 8770/
data (ix(i),i= 401, 480)/
+ 8660, 8580, 8520, 8470, 8440, 8490, 8580, 8660,
1 8710, 8790, 8900, 9070, 8930, 9070, 7840, 7810,
2 7860, 7950, 7950, 8000, 8000, 7920, 7890, 7810,
3 7750, 7700, 7670, 7700, 7670, 7890, 8360, 8410,
4 8440, 8520, 8550, 8600, 8630, 8600, 8630, 8660,
5 8600, 8580, 8580, 8220, 8160, 8110, 8030, 8030,
6 8030, 8080, 8160, 8360, 8220, 8080, 8050, 8140,
7 8140, 8080, 8000, 7920, 7890, 7860, 7860, 7890,
8 7950, 8030, 8140, 8220, 8300, 8190, 8110, 8030,
9 7950, 7920, 7840, 7810, 7810, 7730, 7700, 7700/
data (ix(i),i= 481, 499)/
+ 7700, 7670, 7670, 7700, 7670, 7640, 7590, 7620,
1 7620, 7640, 7670, 7730, 7810, 7950, 8030, 8140,
2 8140, 8220, 8250/
data (iy(i),i= 1, 80)/
+ 10000, 10000, 10000, 0, 0, 0, 90, 120,
1 240, 300, 420, 520, 700, 820, 970, 1060,
2 1090, 1120, 1090, 1120, 1120, 1060, 1060, 1090,
3 1210, 1240, 1240, 1270, 1300, 1270, 1360, 1360,
4 1300, 1270, 1240, 1150, 1000, 880, 850, 820,
5 700, 550, 480, 330, 210, 150, 120, 0,
6 0, 90, 90, 120, 120, 240, 360, 450,
7 580, 670, 700, 640, 580, 480, 390, 330,
8 300, 240, 210, 150, 0, 0, 120, 180,
9 270, 450, 670, 790, 910, 1030, 1180, 1330/
data (iy(i),i= 81, 160)/
+ 1450, 1580, 1760, 1940, 2060, 2150, 2270, 2330,
1 2390, 2390, 2300, 2210, 2120, 1970, 1970, 2000,
2 2030, 2060, 2060, 2060, 2420, 2480, 2580, 2580,
3 2670, 2700, 2760, 2820, 2880, 2970, 3030, 3120,
4 3210, 3270, 3330, 3360, 3330, 3300, 3270, 3300,
5 3330, 3360, 3580, 3640, 3700, 3760, 3790, 3880,
6 3940, 3940, 3910, 3910, 3940, 3970, 4030, 4090,
7 4150, 4180, 4300, 4390, 4450, 4610, 4700, 4880,
8 4910, 4910, 4880, 4880, 4910, 4940, 4940, 4910,
9 4880, 4790, 4790, 4850, 4940, 5060, 5180, 5240/
data (iy(i),i= 161, 240)/
+ 5240, 5270, 5300, 5330, 5390, 5420, 5420, 5520,
1 5580, 5640, 5640, 5670, 5700, 5730, 5670, 5670,
2 5670, 5730, 5880, 5910, 5970, 5970, 5940, 6420,
3 6520, 6610, 6640, 6700, 6850, 6910, 6910, 6910,
4 6970, 7090, 7210, 7240, 7300, 7390, 7580, 7670,
5 7730, 7790, 7880, 7970, 8120, 8150, 8240, 8420,
6 8480, 8520, 8700, 8730, 8790, 8850, 8850, 9150,
7 9150, 9300, 9610, 9700, 9700, 6550, 6420, 6210,
8 6180, 6180, 6300, 6480, 6480, 6580, 6610, 3580,
9 3450, 3330, 3520, 2850, 2760, 2820, 2910, 2640/
data (iy(i),i= 241, 320)/
+ 2580, 2550, 2670, 2090, 2090, 2060, 1970, 2000,
1 2120, 1420, 1360, 1360, 1270, 1270, 1360, 1450,
2 1790, 1730, 1760, 1700, 1790, 1820, 1730, 1700,
3 1580, 1580, 1610, 1580, 1550, 1550, 1610, 1730,
4 1450, 1550, 1580, 1480, 1480, 1520, 1520, 1480,
5 1520, 1390, 1360, 1420, 1330, 1330, 1360, 1450,
6 1390, 1330, 1300, 1270, 1210, 1240, 1240, 1390,
7 1120, 1150, 1150, 1090, 1030, 1000, 1000, 1090,
8 940, 850, 760, 760, 820, 910, 610, 730,
9 730, 640, 610, 580, 550, 480, 480, 390/
data (iy(i),i= 321, 400)/
+ 390, 520, 610, 240, 210, 120, 120, 300,
1 300, 300, 300, 270, 270, 180, 180, 210,
2 180, 240, 300, 150, 90, 90, 60, 60,
3 150, 4270, 4220, 4210, 4180, 4180, 4330, 4270,
4 6390, 6390, 6360, 6210, 6180, 6060, 6000, 6000,
5 6030, 6060, 6150, 6240, 6270, 6300, 6390, 6390,
6 6000, 5970, 6060, 6060, 6000, 5940, 5940, 5880,
7 5820, 5730, 5760, 5880, 6000, 6030, 6450, 6300,
8 6240, 6330, 6480, 6700, 6670, 6640, 6610, 6670,
9 6730, 6640, 6640, 6790, 7210, 7150, 7120, 7090/
data (iy(i),i= 401, 480)/
+ 7090, 7060, 7000, 7030, 7090, 7150, 7150, 7180,
1 7180, 7210, 7240, 7330, 7360, 7420, 7480, 7390,
2 7300, 7150, 7060, 7060, 7000, 7000, 6970, 7030,
3 7090, 7180, 7270, 7300, 7390, 7610, 7820, 7820,
4 7880, 7880, 7850, 7850, 7790, 7730, 7670, 7640,
5 7610, 7520, 7390, 7390, 7480, 7480, 7520, 7640,
6 7730, 7820, 7850, 7850, 9730, 9700, 9640, 9610,
7 9550, 9520, 9450, 9360, 9270, 9210, 9030, 8850,
8 8760, 8700, 8640, 8640, 8640, 8550, 8520, 8450,
9 8520, 8580, 8670, 8640, 8550, 8520, 8480, 8330/
data (iy(i),i= 481, 499)/
+ 8180, 8090, 8030, 8030, 7880, 7760, 7850, 8060,
1 8120, 8700, 8850, 9060, 9330, 9640, 9790, 9790,
2 9760, 9760, 9730/
c
sp(2)='north sea'
sp(1)='north sea'
sp(3)='north sea'
c
do i=1,nvf
vx(i)=real(ix(i),rknd)*1.0e-4_rknd
vy(i)=real(iy(i),rknd)*1.0e-4_rknd
enddo
do i=1,499
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)=0
enddo
ibndry(2,221)=1
ibndry(2,231)=222
ibndry(2,235)=232
ibndry(2,239)=236
ibndry(2,243)=240
ibndry(2,249)=244
ibndry(2,256)=250
ibndry(2,262)=257
ibndry(2,272)=263
ibndry(2,283)=273
ibndry(2,288)=284
ibndry(2,296)=289
ibndry(2,304)=297
ibndry(2,310)=305
ibndry(2,323)=311
ibndry(2,330)=324
ibndry(2,339)=331
ibndry(2,345)=340
ibndry(2,352)=346
ibndry(2,368)=353
ibndry(2,372)=369
ibndry(2,382)=373
ibndry(2,387)=383
ibndry(2,393)=388
ibndry(2,396)=394
ibndry(2,411)=397
ibndry(2,414)=412
ibndry(2,430)=415
ibndry(2,452)=431
ibndry(2,499)=453
do i=500,531
ibndry(3,i)=0
ibndry(4,i)=0
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=0
enddo
ibndry(1,500)=5
ibndry(2,500)=344
ibndry(1,501)=340
ibndry(2,501)=337
ibndry(1,502)=334
ibndry(2,502)=327
ibndry(1,503)=330
ibndry(2,503)=321
ibndry(1,504)=313
ibndry(2,504)=308
ibndry(1,505)=305
ibndry(2,505)=303
ibndry(1,506)=299
ibndry(2,506)=295
ibndry(1,507)=289
ibndry(2,507)=287
ibndry(1,508)=285
ibndry(2,508)=283
ibndry(1,509)=276
ibndry(2,509)=270
ibndry(1,510)=264
ibndry(2,510)=261
ibndry(1,511)=258
ibndry(2,511)=255
ibndry(1,512)=250
ibndry(2,512)=248
ibndry(1,513)=245
ibndry(2,513)=242
ibndry(1,514)=243
ibndry(2,514)=236
ibndry(1,515)=238
ibndry(2,515)=349
ibndry(1,516)=239
ibndry(2,516)=235
ibndry(1,517)=234
ibndry(2,517)=114
ibndry(1,518)=355
ibndry(2,518)=184
ibndry(1,519)=362
ibndry(2,519)=224
ibndry(1,520)=226
ibndry(2,520)=372
ibndry(1,521)=370
ibndry(2,521)=376
ibndry(1,522)=382
ibndry(2,522)=385
ibndry(1,523)=387
ibndry(2,523)=395
ibndry(1,524)=396
ibndry(2,524)=393
ibndry(1,525)=388
ibndry(2,525)=403
ibndry(1,526)=411
ibndry(2,526)=412
ibndry(1,527)=413
ibndry(2,527)=443
ibndry(1,528)=448
ibndry(2,528)=430
ibndry(1,529)=429
ibndry(2,529)=486
ibndry(1,530)=2
ibndry(2,530)=495
ibndry(1,531)=351
ibndry(2,531)=379
c
do i=1,nbf
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd8(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(242) :: ix,iy
integer(kind=iknd), save :: nvf,ntf,nbf
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 domain in the shape of a three piece airfoil
c the outer boundary is an irregular polygon
c data for this domain was provided by tony chan from ucla
c
external sxy
data nvf,ntf,nbf/242,2,246/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
data (ix(i),i= 1, 50)/
+ 1423869, 1423991, 1424665, 1425178, 1425859,
1 1426755, 1427902, 1429365, 1431152, 1433195,
2 1435414, 1437821, 1440377, 1443041, 1445778,
3 1448556, 1451382, 1459973, 1837655, 1898567,
4 1967503, 2440690, 3106079, 3833284, 4532788,
5 5195898, 5774643, 5842507, 5900000, 5951276,
6 5988782, 6013725, 6032087, 6035385, 6029606,
7 6030057, 6035735, 6040665, 6046612, 6053578,
8 6061420, 6070220, 6080025, 6090760, 6102237,
9 6114294, 6126850, 6139864, 6153344, 6167255/
data (ix(i),i= 51, 100)/
+ 6181578, 6196272, 6211289, 6242031, 6257581,
1 6273259, 6289045, 6304912, 6321770, 6392878,
2 6423755, 6417585, 6407490, 6377786, 6144193,
3 5309548, 4401639, 3504295, 2471459, 2220550,
4 2145875, 2075516, 2009864, 1947819, 1889747,
5 1835234, 1784572, 1692365, 1612390, 1577555,
6 1547325, 1521174, 1497263, 1475941, 1457835,
7 1443479, 1433836, 1428603, 1425751, 1424517,
8 1423998, 6435090, 6435684, 6436864, 6438792,
9 6441638, 6445023, 6448381, 6451503, 6453998/
data (ix(i),i= 101, 150)/
+ 6456095, 6500034, 6519509, 6535504, 6544567,
1 6554438, 6565221, 6589578, 6651070, 6689324,
2 6709135, 6727578, 6800035, 6810905, 6844863,
3 6851923, 6858552, 6860805, 6860022, 6858930,
4 6857467, 6855525, 6842807, 6832409, 6818113,
5 6779353, 6725953, 6697470, 6668447, 6636241,
6 6599383, 6562072, 6528875, 6500614, 6477822,
7 6461686, 6452649, 6448563, 6446756, 6445690,
8 6444807, 6438994, 6436979, 6435461, 6613610,
9 6614642, 6619005, 6625937, 6635061, 6690099/
data (ix(i),i= 151, 200)/
+ 6901211, 7152194, 7376001, 7768795, 7782484,
1 7790374, 7788553, 7783294, 7771088, 7137460,
2 7087985, 7039049, 7014894, 6990902, 6967044,
3 6943551, 6920739, 6899039, 6878715, 6859989,
4 6845807, 6832375, 6818686, 6804809, 6790810,
5 6776770, 6762792, 6748974, 6735394, 6722137,
6 6709277, 6696883, 6685029, 6673808, 6663302,
7 6653659, 6644960, 6637259, 6630630, 6625094,
8 6620691, 6617758, 6615033, 0, 439662,
9 679643, 137097, 451695, 1296440, 2369967/
data (ix(i),i= 201, 242)/
+ 3584674, 4993196, 5491886, 6259193, 7002156,
1 8059564, 9048663, 9714154, 9858797, 9717189,
2 9646991, 9421564, 9776057, 9837757, 9778552,
3 9474904, 9602283, 10000000, 9838499, 9353390,
4 8462743, 7620040, 7015043, 6060112, 5380730,
5 4499139, 3494099, 2738492, 1768388, 859205,
6 233243, 251335, 160179, 272935, 94844,
7 248004, 178170, 291858, 115478, 247026,
8 53696, 1598/
data (iy(i),i= 1, 50)/
+ 5254964, 5248606, 5238330, 5233882, 5229771,
1 5225935, 5222397, 5219099, 5216099, 5213480,
2 5211248, 5209293, 5207607, 5206144, 5204861,
3 5203726, 5202719, 5200189, 5137385, 5125430,
4 5115393, 5066633, 5035633, 5050275, 5103050,
5 5192077, 5298739, 5308269, 5311492, 5316925,
6 5327238, 5338626, 5351327, 5377162, 5413147,
7 5443720, 5469391, 5485436, 5501084, 5516263,
8 5530949, 5545001, 5558247, 5570642, 5582239,
9 5593142, 5603388, 5612962, 5621775, 5629776/
data (iy(i),i= 51, 100)/
+ 5636910, 5643118, 5648364, 5656289, 5659461,
1 5661924, 5663409, 5663775, 5663349, 5659469,
2 5659043, 5661356, 5663991, 5670498, 5708333,
3 5810024, 5862396, 5861083, 5800256, 5769966,
4 5754943, 5736496, 5716340, 5694184, 5670140,
5 5644723, 5618731, 5565919, 5513785, 5488443,
6 5464458, 5441265, 5417161, 5391871, 5365872,
7 5339320, 5315813, 5297097, 5282008, 5271090,
8 5262330, 5560329, 5553057, 5545997, 5539125,
9 5533072, 5528495, 5525126, 5522541, 5520937/
data (iy(i),i= 101, 150)/
+ 5519998, 5501572, 5494669, 5489811, 5487370,
1 5485138, 5483336, 5480422, 5469531, 5460928,
2 5455621, 5449801, 5423853, 5418400, 5395982,
3 5391642, 5387869, 5389659, 5392076, 5394584,
4 5397235, 5400257, 5417443, 5430332, 5446495,
5 5486280, 5535828, 5559779, 5581622, 5601534,
6 5618159, 5628927, 5632761, 5630887, 5625092,
7 5617436, 5610627, 5606162, 5603417, 5601282,
8 5599234, 5582132, 5575334, 5567821, 5167153,
9 5156108, 5145026, 5134288, 5124168, 5077804/
data (iy(i),i= 151, 200)/
+ 4898398, 4659643, 4419134, 3942617, 3929639,
1 3923042, 3933236, 3951430, 3975541, 4989411,
2 5054514, 5112266, 5138064, 5161310, 5181898,
3 5199488, 5214413, 5226979, 5237375, 5245711,
4 5251225, 5255710, 5259456, 5262411, 5264493,
5 5265598, 5265665, 5264710, 5262919, 5260344,
6 5257038, 5253061, 5248457, 5243203, 5237352,
7 5230837, 5223771, 5216286, 5208517, 5200622,
8 5192757, 5186040, 5177227, 2881950, 2592188,
9 1793022, 1106241, 258165, 210040, 321775/
data (iy(i),i= 201, 242)/
+ 0, 413863, 443966, 656610, 223114,
1 139355, 589955, 1421292, 2288002, 3260701,
2 3797993, 4254527, 4765088, 5251086, 5694451,
3 6190597, 6807999, 7503138, 8281126, 9070889,
4 9623615, 9899876, 9228654, 9173282, 9671695,
5 10000000, 9896622, 9445850, 9858504, 9651270,
6 9304286, 8609588, 7904434, 7103745, 6408539,
7 5841843, 5419985, 5018809, 4774224, 4486439,
8 4057561, 3433861/
c
sp(2)='airfoil'
sp(1)='airfoil'
sp(3)='airfoil'
c
do i=1,nvf
vx(i)=real(ix(i),rknd)*1.0e-7_rknd
vy(i)=real(iy(i),rknd)*1.0e-7_rknd
enddo
do i=1,nbf
ibndry(1,i)=i
ibndry(2,i)=i+1
ibndry(3,i)=0
ibndry(4,i)=2
if(i>nvf) ibndry(4,i)=0
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=1
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
ibndry(2,91)=1
ibndry(2,144)=92
ibndry(2,193)=145
ibndry(2,242)=194
c
ibndry(1,243)=237
ibndry(2,243)=2
ibndry(1,244)=61
ibndry(2,244)=136
ibndry(1,245)=117
ibndry(2,245)=172
ibndry(1,246)=156
ibndry(2,246)=212
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd9(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(39) :: ib1,ib2,ib3,ib4
integer(kind=iknd), save, dimension(30) :: ix,iy
integer(kind=iknd), save :: nvf,ntf,nbf
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(100) :: rp,ru
character(len=80), dimension(100) :: sp,su
cy
c a domain composed of nine subregions
c data for this problem was supplied by hans mittelmann
c
external sxy
data nvf,ntf,nbf/30,9,39/
data (ix(i),i= 1, 30)/
+ -3000,-3000,-3000,-3000, 3600, 7600, 4100, 7100, 4700,
1 7100, 3600, 4100, 4700, 7100, 9000, 9000, 9000, 7600,
2 575, 775, 0, 0, -575, -775, 0, 0, 500,
3 0, -500, 0/
data (iy(i),i= 1, 30)/
+ 0,-4400, 3200, 3600, 0, 0, 500, 500, 1100,
1 1100, 3200, 3600, 3600, 3600, 3600, 3200,-4400, 3200,
2 0, 0, 575, 775, 0, 0, -575, -775, 0,
3 500, 0, -500/
data ib1/ 4,12,13,14, 3,18, 9, 7, 5, 4, 3,11,12,
+ 13,14,18,15,16,10, 2, 1, 1,24,24,22,20,
1 20,26,23,21,19,25,19,23,29,27,27,28,19/
data ib2/12,13,14,15,11,16,10, 8, 6, 3, 1, 5, 7,
+ 9,10, 6,16,17, 8,17, 2,24,23,22,20, 5,
1 26,24,21,19,25,23,20,29,30,30,28,29,27/
data ib3/0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
+ 1,1,0,1,1,1,1,1,1,0,0,1,1,1,1,0/
data ib4/1,1,1,1,0,0,0,0,0,1,1,0,0,0,0,0,1,1,0,1,1,0,0,0,0,
+ 0,0,0,0,0,0,0,0,0,2,2,2,2,0/
c
sp(2)='plant'
sp(1)='plant'
sp(3)='plant'
c
do i=1,nbf
ibndry(1,i)=ib1(i)
ibndry(2,i)=ib2(i)
ibndry(3,i)=ib3(i)
ibndry(4,i)=ib4(i)
ibndry(5,i)=0
ibndry(6,i)=0
if(i<=4) then
ibndry(7,i)=1
else
ibndry(7,i)=0
endif
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
do i=1,nvf
vx(i)=real(ix(i),rknd)*1.0e-2_rknd
vy(i)=real(iy(i),rknd)*1.0e-2_rknd
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
c
c make itnode, divide long curved edges
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
call sklutl(1_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
c
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd10(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
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(100) :: rp,ru,xm,ym
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c circular region with many curved interior edges
c dirichlet boundary conditions on all boundary edges
c
external sxy
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
sp(2)='fan'
sp(1)='fan'
sp(3)='fan'
c
pi=3.141592653589793e0_rknd
r=1.0e0_rknd
r2=r/sqrt(2.0e0_rknd)
ntf=7
nvf=ntf+1
nbf=2*ntf
vx(1)=0.0e0_rknd
vy(1)=0.0e0_rknd
dt=2.0e0_rknd*pi/real(ntf,rknd)
xm(1)=0.0e0_rknd
ym(1)=0.0e0_rknd
do i=1,ntf
ang=real(i-1,rknd)*dt
angm=ang-pi/4.0e0_rknd
vx(i+1)=r*cos(ang)
vy(i+1)=r*sin(ang)
xm(i+1)=r2*cos(angm)
ym(i+1)=r2*sin(angm)
k=2*i-1
ibndry(1,k)=1
ibndry(2,k)=i+1
ibndry(3,k)=i+1
ibndry(4,k)=0
ibndry(5,k)=0
ibndry(6,k)=0
ibndry(7,k)=k
ibndry(1,k+1)=i+1
ibndry(2,k+1)=i+2
if(i==ntf) ibndry(2,k+1)=2
ibndry(3,k+1)=1
ibndry(4,k+1)=2
ibndry(5,k+1)=0
ibndry(6,k+1)=0
ibndry(7,k+1)=k+1
enddo
c
do i=1,nbf
if(ibndry(3,i)>0) then
sf(1,i)=xm(ibndry(3,i))
sf(2,i)=ym(ibndry(3,i))
else
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
endif
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode, refine long edges, find symmetries
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
call sklutl(1_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
call sklutl(2_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd11(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
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c domain which is decomposed into several distinct pieces,
c and making extensive use of linked edges.
c
external sxy
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
sp(2)='arc'
sp(1)='arc'
sp(3)='arc'
c
pi=3.141592653589793e0_rknd
r1=1.0e0_rknd
r2=2.0e0_rknd
ntf=6
nvf=4*ntf
nbf=nvf
dt=pi/real(ntf,rknd)
eps=dt/20.0e0_rknd
c
do i=1,nbf
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
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
do i=1,ntf
k=(i-1)*4
t1=real(i-1,rknd)*dt+eps
t2=real(i,rknd)*dt-eps
c1=cos(t1)
c2=cos(t2)
s1=sin(t1)
s2=sin(t2)
vx(k+1)=r1*c1
vy(k+1)=r1*s1
vx(k+2)=r2*c1
vy(k+2)=r2*s1
vx(k+3)=r2*c2
vy(k+3)=r2*s2
vx(k+4)=r1*c2
vy(k+4)=r1*s2
ibndry(2,k+4)=k+1
ibndry(3,k+2)=1
ibndry(3,k+4)=1
if(i>1) ibndry(4,k+1)=-(k-1)
if(i<ntf) ibndry(4,k+3)=-(k+5)
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode, refine long edges, find symmetries
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
call sklutl(1_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
call sklutl(2_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd12(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 :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru,xm,ym
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c a spiral shaped domain for testing tgen
c
external sxy
data ntf,nvf,nbf/1,0,0/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
sp(2)='spiral'
sp(1)='spiral'
sp(3)='spiral'
c
nloops=6
nvf=8*nloops+2
nbf=nvf
c
xm(1)=0.0e0_rknd
ym(1)=0.0e0_rknd
xm(2)=1.0e0_rknd
ym(2)=0.0e0_rknd
c
do k=1,nbf
ibndry(1,k)=k
ibndry(2,k)=k+2
ibndry(3,k)=1
ibndry(4,k)=2
ibndry(5,k)=0
ibndry(6,k)=0
ibndry(7,k)=k
enddo
c
do i=1,nloops
k=(i-1)*8
m=2*i+1
vx(k+1)=-real(m-1,rknd)
vy(k+1)=0.0e0_rknd
vx(k+2)=-real(m,rknd)
vy(k+2)=0.0e0_rknd
vx(k+3)=1.0e0_rknd
vy(k+3)=real(m,rknd)
vx(k+4)=1.0e0_rknd
vy(k+4)=real(m+1,rknd)
vx(k+5)=real(m+1,rknd)
vy(k+5)=0.0e0_rknd
vx(k+6)=real(m+2,rknd)
vy(k+6)=0.0e0_rknd
vx(k+7)=0.0e0_rknd
vy(k+7)=-real(m+1,rknd)
vx(k+8)=0.0e0_rknd
vy(k+8)=-real(m+2,rknd)
do j=k+1,k+4
ibndry(3,j)=2
enddo
enddo
m=2*nloops+3
k=nloops*8
vx(k+1)=-real(m-1,rknd)
vy(k+1)=0.0e0_rknd
vx(k+2)=-real(m,rknd)
vy(k+2)=0.0e0_rknd
ibndry(1,k+1)=1
ibndry(2,k+1)=2
ibndry(3,k+1)=0
ibndry(1,k+2)=k+1
ibndry(2,k+2)=k+2
ibndry(3,k+2)=0
c
itnode(1,1)=1
itnode(2,1)=1
itnode(3,1)=0
itnode(5,1)=1
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
do i=1,nbf
if(ibndry(3,i)>0) then
sf(1,i)=xm(ibndry(3,i))
sf(2,i)=ym(ibndry(3,i))
else
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
endif
enddo
c
c refine long edges
c
call sklutl(1_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
c
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd13(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(135) :: ix,iy
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c this domain is a facsimile of the newly adopted ucsd logo.
c dirichlet boundary conditions are imposed.
c
external sxy
data (ix(i),i= 1, 90)/
+ 320, 520, 520, 450, 460, 320, 250, 585, 585,
1 460, 470, 250, 165, 658, 658, 470, 480, 165,
2 100, 680, 680, 600, 550, 520, 490, 470, 460,
3 100, 460, 500, 530, 550, 563, 570, 560, 540,
4 510, 280, 310, 316, 320, 319, 316, 310, 100,
5 40, 68, 68, 75, 100, 200, 220, 235, 235,
6 190, 215, 215, 205, 190, 120, 105, 100, 390,
7 390, 370, 310, 280, 255, 245, 245, 255, 270,
8 300, 380, 390, 330, 300, 285, 270, 270, 280,
9 295, 310, 350, 380, 545, 545, 525, 445, 420/
data (ix(i),i= 91, 135)/
+ 410, 410, 420, 500, 520, 500, 420, 420, 510,
1 540, 555, 555, 545, 500, 445, 435, 450, 505,
2 535, 610, 550, 578, 578, 610, 610, 670, 700,
3 715, 715, 700, 670, 620, 620, 690, 725, 737,
4 747, 747, 737, 725, 690, -332, 1147, 1147, -332/
data (iy(i),i= 1, 90)/
+ 520, 520, 490, 490, 500, 500, 482, 482, 450,
1 450, 460, 460, 440, 440, 410, 410, 420, 420,
2 400, 400, 370, 360, 350, 345, 330, 310, 290,
3 380, 270, 260, 250, 240, 230, 219, 210, 207,
4 210, 280, 281, 283, 288, 293, 296, 300, 150,
5 150, 140, 30, 10, 0, 0, 13, 25, 150,
6 150, 140, 30, 20, 10, 10, 20, 30, 110,
7 145, 150, 150, 145, 130, 100, 50, 30, 15,
8 0, 0, 10, 10, 20, 30, 50, 100, 120,
9 130, 140, 140, 115, 110, 145, 150, 150, 140/
data (iy(i),i= 91, 135)/
+ 130, 100, 85, 50, 30, 10, 10, 0, 0,
1 10, 25, 50, 60, 80, 100, 120, 140, 140,
2 115, 150, 150, 140, 0, 0, 140, 140, 130,
3 105, 45, 20, 10, 10, 0, 0, 15, 30,
4 50, 100, 120, 135, 150, -400, -400, 920, 920/
c
data ntf,nvf,nbf/9,135,145/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
sp(2)='ucsd'
sp(1)='ucsd'
sp(3)='ucsd'
c
do i=1,nvf
vx(i)=real(ix(i),rknd)
vy(i)=real(iy(i),rknd)
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
enddo
ibndry(2,6)=1
ibndry(2,12)=7
ibndry(2,18)=13
ibndry(2,27)=19
ibndry(2,44)=28
ibndry(2,62)=45
ibndry(2,85)=63
ibndry(2,109)=86
ibndry(2,131)=110
ibndry(2,135)=132
c
do i=nvf+1,nbf
ibndry(3,i)=0
ibndry(4,i)=0
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=i
enddo
ibndry(1,136)=132
ibndry(2,136)=49
ibndry(1,137)=54
ibndry(2,137)=68
ibndry(1,138)=75
ibndry(2,138)=97
ibndry(1,139)=100
ibndry(2,139)=113
ibndry(1,140)=110
ibndry(2,140)=35
ibndry(1,141)=29
ibndry(2,141)=27
ibndry(1,142)=19
ibndry(2,142)=18
ibndry(1,143)=13
ibndry(2,143)=12
ibndry(1,144)=6
ibndry(2,144)=7
ibndry(1,145)=1
ibndry(2,145)=135
c
do i=1,nbf
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd14(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 :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(16) :: xx,yy
character(len=80), dimension(100) :: sp,su
cy
external sxy
data ntf,nvf,nbf/4,29,32/
c
c a domain in the shape of a nozzle.
c four similar regions, several curved edges.
c
sp(2)='nozzle'
sp(1)='nozzle'
sp(3)='nozzle'
c
rp(1)=ru(8)
c
vx(1)=0.0e0_rknd
vy(1)=0.0e0_rknd
vx(2)=4.0e0_rknd
vy(2)=0.0e0_rknd
vx(3)=4.0e0_rknd
vy(3)=1.0e0_rknd
vx(4)=3.0e0_rknd
vy(4)=1.0e0_rknd
vx(5)=2.5e0_rknd
vy(5)=f(2.5e0_rknd)
vx(6)=2.0e0_rknd
vy(6)=f(2.0e0_rknd)
vx(7)=1.5e0_rknd
vy(7)=f(1.5e0_rknd)
vx(8)=1.0e0_rknd
vy(8)=0.5e0_rknd
vx(9)=0.0e0_rknd
vy(9)=0.5e0_rknd
do i=10,16
vx(i)=-vx(18-i)
vy(i)=vy(18-i)
enddo
do i=17,29
vx(i)=vx(32-i)
vy(i)=-vy(32-i)
enddo
c
xx(1)=2.75e0_rknd
yy(1)=f(2.75e0_rknd)
xx(2)=2.25e0_rknd
yy(2)=f(2.25e0_rknd)
xx(3)=1.75e0_rknd
yy(3)=f(1.75e0_rknd)
xx(4)=1.25e0_rknd
yy(4)=f(1.25e0_rknd)
do i=1,4
xx(4+i)=-xx(5-i)
yy(4+i)=yy(5-i)
xx(8+i)=-xx(i)
yy(8+i)=-yy(i)
xx(12+i)=xx(5-i)
yy(12+i)=-yy(5-i)
enddo
c
do i=1,29
ibndry(1,i)=i
ibndry(2,i)=i+1
ibndry(3,i)=0
ibndry(4,i)=1
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=i
enddo
do i=1,4
ibndry(3,3+i)=i
ibndry(3,9+i)=i+4
ibndry(3,17+i)=i+8
ibndry(3,23+i)=i+12
enddo
ibndry(4,1)=0
ibndry(4,2)=2
ibndry(4,15)=2
ibndry(4,16)=2
ibndry(4,29)=2
ibndry(2,29)=2
do i=30,32
ibndry(1,i)=1
ibndry(2,i)=9
ibndry(3,i)=0
ibndry(4,i)=0
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=i
enddo
ibndry(2,30)=16
ibndry(2,32)=23
c
do i=1,nbf
if(ibndry(3,i)>0) then
i1=ibndry(1,i)
i2=ibndry(2,i)
i3=ibndry(3,i)
call centre(vx(i1),vy(i1),vx(i2),vy(i2),
+ xx(i3),yy(i3),sf(1,i),sf(2,i))
else
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
endif
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
c
c make itnode, find symmetries
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
call sklutl(2_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
c
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
function f(x)
cx
use mthdef
implicit real(kind=rknd) (a-h,o-z)
implicit integer(kind=iknd) (i-n)
real(kind=rknd) :: f
cy
f=-0.25e0_rknd*x+1.5e0_rknd
f=f*x-2.25e0_rknd
f=(f*x+2.0e0_rknd)/2.0e0_rknd
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd15(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(170) :: ix,iy
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
external sxy
data (ix(i),i= 1, 90)/
+ 0, 40, 80, 110, 130, 140, 90, 60, 100,
1 180, 230, 210, 260, 300, 350, 380, 370, 330,
2 400, 430, 425, 420, 430, 455, 490, 525, 560,
3 580, 590, 585, 570, 540, 500, 530, 570, 610,
4 630, 580, 610, 650, 690, 730, 715, 775, 755,
5 700, 700, 710, 740, 760, 790, 815, 830, 890,
6 865, 900, 950, 975, 1000, 1040, 1075, 1085, 1100,
7 1130, 1115, 1115, 1130, 1170, 1200, 1230, 1250, 1250,
8 1300, 1310, 1300, 1290, 1420, 1360, 1270, 1240, 1170,
9 1100, 1000, 900, 800, 700, 600, 500, 400, 350/
data (ix(i),i= 91, 170)/
+ 320, 310, 320, 350, 400, 500, 600, 700, 800,
1 900, 1000, 1100, 1200, 1160, 1100, 1050, 1025, 1020,
2 1010, 1000, 980, 950, 895, 905, 915, 905, 875,
3 840, 780, 790, 750, 700, 660, 620, 570, 520,
4 470, 440, 400, 350, 300, 240, 180, 150, 110,
5 60, 190, 220, 280, 310, 280, 220, 480, 485,
6 520, 550, 550, 520, 1170, 1200, 1240, 1260, 1220,
7 1175, 1180, 1100, 1000, 900, 800, 700, 600, 500,
8 460, 500, 600, 700, 800, 900, 1000, 1100/
data (iy(i),i= 1, 90)/
+ 70, 100, 160, 200, 260, 300, 270, 310, 340,
1 380, 380, 340, 355, 370, 360, 330, 280, 240,
2 270, 300, 340, 370, 400, 450, 500, 545, 575,
3 550, 505, 450, 400, 350, 300, 265, 280, 320,
4 350, 350, 390, 400, 480, 500, 410, 420, 360,
5 350, 300, 260, 250, 270, 300, 340, 390, 400,
6 360, 385, 400, 380, 390, 390, 370, 300, 270,
7 260, 300, 320, 360, 380, 385, 380, 365, 380,
8 390, 300, 240, 200, 220, 150, 130, 100, 50,
9 10, -5, -15, -10, -5, 0, 15, 45, 70/
data (iy(i),i= 91, 170)/
+ 100, 120, 140, 150, 160, 161, 162, 163, 165,
1 166, 167, 168, 170, 215, 200, 200, 220, 270,
2 315, 330, 300, 210, 200, 240, 290, 310, 290,
3 210, 200, 220, 210, 200, 220, 245, 210, 200,
4 220, 250, 225, 210, 200, 190, 200, 150, 100,
5 60, 260, 290, 320, 310, 280, 250, 350, 400,
6 460, 500, 460, 400, 330, 340, 310, 270, 260,
7 300, 130, 70, 50, 45, 50, 60, 70, 80,
8 95, 110, 120, 122, 124, 126, 127, 129/
c
data ntf,nvf,nbf/5,170,178/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
sp(2)='pltmg'
sp(1)='pltmg'
sp(3)='pltmg'
c
do i=1,nvf
vx(i)=real(ix(i),rknd)
vy(i)=real(iy(i),rknd)
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
enddo
ibndry(2,136)=1
ibndry(2,142)=137
ibndry(2,148)=143
ibndry(2,154)=149
ibndry(2,170)=155
do i=nvf+1,nbf
ibndry(3,i)=0
ibndry(4,i)=0
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=i
enddo
ibndry(1,171)=5
ibndry(2,171)=137
ibndry(1,172)=16
ibndry(2,172)=140
ibndry(1,173)=20
ibndry(2,173)=143
ibndry(1,174)=27
ibndry(2,174)=146
ibndry(1,175)=67
ibndry(2,175)=149
ibndry(1,176)=75
ibndry(2,176)=152
ibndry(1,177)=102
ibndry(2,177)=170
ibndry(1,178)=82
ibndry(2,178)=156
c
do i=1,nbf
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd16(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
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c domain which is has many cracks
c
external sxy
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
sp(2)='crack'
sp(1)='crack'
sp(3)='crack'
c
pi=3.141592653589793e0_rknd
r1=1.0e0_rknd
r2=2.0e0_rknd
ncrack=8
nvf=1+3*ncrack
nbf=4*ncrack
ntf=ncrack
c
vx(1)=0.0e0_rknd
vy(1)=0.0e0_rknd
do i=1,ncrack
ang=2.0e0_rknd*pi*real(i-1,rknd)/real(ncrack,rknd)
c=cos(ang)
s=sin(ang)
ii=3*i
vx(ii-1)=r2*c
vy(ii-1)=r2*s
vx(ii)=r1*c
vy(ii)=r1*s
vx(ii+1)=r2*c
vy(ii+1)=r2*s
c
jj=4*i-3
c
ibndry(1,jj)=1
ibndry(2,jj)=ii
ibndry(3,jj)=0
ibndry(4,jj)=0
ibndry(5,jj)=0
ibndry(6,jj)=0
ibndry(7,jj)=0
c
ibndry(1,jj+1)=ii
ibndry(2,jj+1)=ii-1
ibndry(3,jj+1)=0
ibndry(4,jj+1)=1
ibndry(5,jj+1)=0
ibndry(6,jj+1)=0
ibndry(7,jj+1)=1
c
ibndry(1,jj+2)=ii
ibndry(2,jj+2)=ii+1
ibndry(3,jj+2)=0
ibndry(4,jj+2)=2
ibndry(5,jj+2)=0
ibndry(6,jj+2)=0
ibndry(7,jj+2)=2
c
ibndry(1,jj+3)=ii+1
ibndry(2,jj+3)=ii+2
ibndry(3,jj+3)=1
ibndry(4,jj+3)=2
ibndry(5,jj+3)=0
ibndry(6,jj+3)=0
ibndry(7,jj+3)=3
enddo
ibndry(2,nbf)=2
c
do i=1,nbf
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode, refine long arcs, find symmetries
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
call sklutl(1_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
call sklutl(2_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd17(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), dimension(1250) :: list
integer(kind=iknd), save, dimension(300) :: ix,iy
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c region in the shape of monterey bay
c the data was provided by francios lekien
c
external sxy
data ntf,nvf,nbf/1,300,300/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
data (ix(i),i= 1, 50)/
+ -122324501, -122324501, -122135002, -122135002, -122134003,
1 -122133003, -122128998, -122127998, -122125999, -122125000,
2 -122122002, -122120003, -122120003, -122116997, -122115997,
3 -122113998, -122113998, -122111000, -122109001, -122107002,
4 -122105003, -122102997, -122100998, -122098999, -122097000,
5 -122094002, -122089996, -122088997, -122086998, -122084999,
6 -122083000, -122081001, -122080002, -122078003, -122075996,
7 -122075996, -122073997, -122070999, -122069000, -122065002,
8 -122060997, -122059998, -122057999, -122056000, -122052002,
9 -122050003, -122046997, -122044998, -122042000, -122040001/
data (ix(i),i= 51, 100)/
+ -122039001, -122035004, -122033997, -122030998, -122028000,
1 -122024002, -122024002, -122026001, -122026001, -122025002,
2 -122022003, -122016998, -122014000, -122012001, -122008003,
3 -122004997, -122001999, -122001999, -122000000, -121997002,
4 -121992996, -121987999, -121985001, -121980003, -121977997,
5 -121975998, -121973999, -121973000, -121970001, -121968002,
6 -121964996, -121961998, -121959999, -121958000, -121955002,
7 -121953003, -121947998, -121945000, -121941002, -121939003,
8 -121936996, -121933998, -121931999, -121929001, -121925003,
9 -121919998, -121914001, -121911003, -121910004, -121906998/
data (ix(i),i= 101, 150)/
+ -121903000, -121900002, -121897003, -121892998, -121891998,
1 -121889000, -121887001, -121884003, -121883003, -121882004,
2 -121879997, -121876999, -121876999, -121875000, -121872002,
3 -121869003, -121865997, -121862999, -121859001, -121857002,
4 -121855003, -121852997, -121849998, -121848000, -121847000,
5 -121843002, -121838997, -121835999, -121834000, -121830002,
6 -121828003, -121825996, -121822998, -121820999, -121819000,
7 -121816002, -121814003, -121811996, -121808998, -121807999,
8 -121807999, -121806000, -121805000, -121804001, -121804001,
9 -121803001, -121800003, -121799004, -121799004, -121796997/
data (ix(i),i= 151, 200)/
+ -121793999, -121792000, -121790001, -121790001, -121789001,
1 -121789001, -121788002, -121788002, -121789001, -121791000,
2 -121792999, -121793999, -121794998, -121794998, -121795998,
3 -121796997, -121799004, -121800003, -121800003, -121801003,
4 -121802002, -121804001, -121805000, -121805000, -121804001,
5 -121805000, -121806999, -121806999, -121806999, -121806999,
6 -121806999, -121807999, -121808998, -121810997, -121813004,
7 -121815002, -121817001, -121819000, -121821999, -121824997,
8 -121829002, -121832001, -121834999, -121838997, -121842003,
9 -121847000, -121851997, -121859001, -121862999, -121867996/
data (ix(i),i= 201, 250)/
+ -121873001, -121876999, -121882004, -121886002, -121889000,
1 -121892998, -121893997, -121893997, -121896004, -121897003,
2 -121898003, -121900002, -121902000, -121903999, -121903999,
3 -121904999, -121906998, -121908997, -121910004, -121910004,
4 -121911003, -121913002, -121914001, -121915001, -121916000,
5 -121915001, -121915001, -121916000, -121917999, -121919998,
6 -121919998, -121921997, -121921997, -121922997, -121924004,
7 -121926003, -121928001, -121930000, -121931000, -121931999,
8 -121934998, -121936996, -121936996, -121935997, -121936996,
9 -121936996, -121939003, -121939003, -121939003, -121941002/
data (ix(i),i= 251, 300)/
+ -121941002, -121942001, -121941002, -121941002, -121942001,
1 -121944000, -121947998, -121950996, -121950996, -121953003,
2 -121954002, -121956001, -121958000, -121958000, -121959000,
3 -121959999, -121959999, -121960999, -121962997, -121960999,
4 -121961998, -121961998, -121963997, -121964996, -121963997,
5 -121963997, -121966003, -121968002, -121970001, -121973999,
6 -121974998, -121974998, -121974998, -121975998, -121975998,
7 -121977997, -121977997, -121975998, -121974998, -121973000,
8 -121972000, -121971001, -121970001, -121969002, -121967003,
9 -121963997, -121962997, -121959999, -121958000, -121957001/
data (iy(i),i= 1, 50)/
+ 36565899, 36973900, 36973900, 36973202, 36971298,
1 36972500, 36972099, 36971100, 36970001, 36970699,
2 36970699, 36969799, 36967499, 36967499, 36967899,
3 36967201, 36965900, 36965500, 36963902, 36962502,
4 36961399, 36961399, 36960701, 36960499, 36960499,
5 36959801, 36959099, 36959099, 36957901, 36958401,
6 36957298, 36957901, 36958599, 36957901, 36957500,
7 36956402, 36956402, 36954700, 36953800, 36953800,
8 36953899, 36954700, 36954498, 36953400, 36954498,
9 36954498, 36954498, 36955200, 36955200, 36957100/
data (iy(i),i= 51, 100)/
+ 36957901, 36957298, 36956402, 36956799, 36956402,
1 36955700, 36958801, 36961300, 36963902, 36966599,
2 36967899, 36967899, 36968601, 36968601, 36967999,
3 36967999, 36966801, 36965199, 36966099, 36966099,
4 36964500, 36962502, 36962502, 36962502, 36960701,
5 36958801, 36958801, 36958801, 36960499, 36962700,
6 36964600, 36965500, 36967999, 36970501, 36973202,
7 36975201, 36976601, 36978001, 36980202, 36981998,
8 36981998, 36981998, 36982700, 36982201, 36980499,
9 36979099, 36977299, 36975399, 36973900, 36973202/
data (iy(i),i= 101, 150)/
+ 36972099, 36970699, 36969101, 36967300, 36965500,
1 36963200, 36961102, 36958801, 36957500, 36957298,
2 36955700, 36953400, 36952000, 36949100, 36947300,
3 36944099, 36941101, 36936100, 36932201, 36928799,
4 36924801, 36921398, 36917999, 36914799, 36911800,
5 36907700, 36902699, 36897999, 36893799, 36889099,
6 36884499, 36880600, 36875900, 36872299, 36868599,
7 36863899, 36860001, 36856998, 36855400, 36852299,
8 36850700, 36848400, 36845001, 36842499, 36839901,
9 36837299, 36835400, 36833099, 36829800, 36826500/
data (iy(i),i= 151, 200)/
+ 36823200, 36820400, 36817299, 36814701, 36811298,
1 36809898, 36807400, 36804501, 36800900, 36797699,
2 36793598, 36790401, 36786701, 36785000, 36781101,
3 36777199, 36772400, 36769001, 36765202, 36761799,
4 36760201, 36757401, 36753601, 36749699, 36748402,
5 36745899, 36741600, 36737202, 36732899, 36727200,
6 36718300, 36709499, 36701302, 36695900, 36689201,
7 36683601, 36676998, 36670601, 36665600, 36658798,
8 36651901, 36648300, 36643600, 36637199, 36632198,
9 36626099, 36620800, 36614300, 36611500, 36609001/
data (iy(i),i= 201, 250)/
+ 36607201, 36606098, 36604900, 36604900, 36605400,
1 36605801, 36606998, 36611000, 36612900, 36615799,
2 36617699, 36619499, 36621498, 36623798, 36624500,
3 36624199, 36624500, 36623600, 36623600, 36624901,
4 36626099, 36626099, 36626999, 36627800, 36628502,
5 36629200, 36629700, 36629700, 36630402, 36631500,
6 36632900, 36633801, 36635799, 36636799, 36637901,
7 36638802, 36639900, 36640400, 36640400, 36640400,
8 36640400, 36637600, 36635101, 36632198, 36631001,
9 36629902, 36629902, 36628101, 36626900, 36626900/
data (iy(i),i= 251, 300)/
+ 36624500, 36623100, 36621700, 36621101, 36620399,
1 36617699, 36614300, 36612900, 36611900, 36611900,
2 36612598, 36614201, 36612598, 36611099, 36609501,
3 36609200, 36607601, 36605099, 36603100, 36600800,
4 36599201, 36597198, 36596100, 36595600, 36593498,
5 36590099, 36588100, 36585800, 36585800, 36585800,
6 36584900, 36583599, 36582401, 36582401, 36583099,
7 36583099, 36581902, 36580799, 36578999, 36578098,
8 36576099, 36574402, 36573101, 36571899, 36571201,
9 36570301, 36568501, 36567799, 36567402, 36565899/
c
c
sp(2)='monterey bay'
sp(1)='monterey bay'
sp(3)='monterey bay'
c
do i=1,nvf
ibndry(1,i)=i-1
ibndry(2,i)=i
ibndry(3,i)=0
ibndry(4,i)=2
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=0
vx(i)=real(ix(i),rknd)*1.0e-6_rknd
vy(i)=real(iy(i),rknd)*1.0e-6_rknd
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
ibndry(1,1)=nvf
c
c make itnode
c
itnode(1,1)=1
itnode(2,1)=1
itnode(3,1)=0
itnode(4,1)=0
itnode(5,1)=1
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd18(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(89) :: ix,iy
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c siam logo
c
external sxy
data ntf,nvf,nbf/2,89,92/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
data (ix(i),i= 1, 89)/
+ 0, 30, 32, 34, 57, 59, 61, 89, 89,
1 91, 107, 109, 109, 114, 116, 116, 125, 125,
2 130, 132, 132, 134, 158, 158, 141, 139, 139,
3 137, 127, 125, 123, 111, 109, 109, 102, 102,
4 100, 91, 89, 89, 87, 61, 61, 82, 82,
5 81, 67, 66, 66, 67, 80, 80, 60, 58,
6 58, 56, 50, 48, 48, 41, 41, 39, 35,
7 33, 33, 31, 9, 8, 8, 9, 33, 33,
8 3, 1, 1, 3, 25, 26, 26, 25, 0,
9 48, 50, 39, 41, -60, 220, 220, -60/
data (iy(i),i= 1, 89)/
+ 0, 0, 2, 0, 0, 2, 0, 0, 2,
1 0, 0, 2, 26, 26, 24, 0, 0, 26,
2 26, 24, 2, 0, 0, 7, 7, 9, 31,
3 33, 33, 31, 33, 33, 31, 33, 33, 9,
4 7, 7, 9, 31, 33, 33, 27, 27, 8,
5 7, 7, 8, 13, 14, 14, 21, 21, 19,
6 9, 7, 7, 9, 33, 33, 9, 7, 7,
7 9, 19, 21, 21, 22, 26, 27, 27, 33,
8 33, 31, 16, 14, 14, 13, 8, 7, 7,
9 40, 51, 51, 40, -90, -90, 140, 140/
c
sp(2)='siam'
sp(1)='siam'
sp(3)='siam'
c
do i=1,nvf
vx(i)=real(ix(i),rknd)
vy(i)=real(iy(i),rknd)
enddo
do i=1,nbf
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)=0
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
ibndry(2,81)=1
ibndry(2,85)=82
ibndry(2,89)=86
ibndry(1,90)=86
ibndry(2,90)=1
ibndry(1,91)=60
ibndry(2,91)=85
ibndry(1,92)=84
ibndry(2,92)=89
do i=nvf+1,nbf
ibndry(4,i)=0
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd19(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), dimension(1000) :: b1,b2
integer(kind=iknd), dimension(968) :: ix,iy
integer(kind=iknd), save :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
c domain in the shape of mexico
c
external sxy
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
c
data ntf,nvf,nbf/1,968,999/
data nbound/ 427/
c
data (ix(i),i= 1, 80)/
+ 10555, 10555, 10699, 11121, 11236, 11562, 11677, 11677,
1 11792, 11869, 12032, 12032, 11955, 12003, 12358, 12378,
2 12493, 12445, 12608, 12646, 12531, 12838, 13088, 13155,
3 13472, 13587, 13481, 13596, 14057, 14172, 14191, 14373,
4 14776, 14959, 15026, 15170, 15573, 15410, 15007, 14796,
5 14661, 14383, 14162, 13855, 13529, 12579, 11860, 11764,
6 11716, 11600, 11437, 11303, 11140, 11351, 10948, 10699,
7 9355, 9211, 8962, 8616, 8405, 8146, 7686, 7590,
8 6870, 6937, 7446, 7858, 8041, 8252, 8530, 8578,
9 8501, 8300, 8070, 8271, 8041, 7954, 8261, 13951/
data (ix(i),i= 81, 160)/
+ 17396, 16168, 16129, 15832, 15582, 15285, 15131, 14767,
1 14508, 14182, 14066, 14134, 13731, 13644, 13596, 13337,
2 13040, 12790, 12224, 12042, 12042, 11668, 11620, 11389,
3 11399, 10862, 10842, 10670, 10411, 10497, 10497, 10200,
4 9989, 9451, 9221, 8943, 8396, 26300, 18307, 18394,
5 18413, 21138, 21570, 21944, 22548, 22961, 23642, 23825,
6 24016, 24045, 24419, 24698, 25062, 25341, 25696, 25907,
7 32057, 31635, 31280, 31213, 30599, 29975, 29226, 28651,
8 27960, 27749, 27519, 27231, 27020, 26549, 4299, 8146,
9 8002, 7839, 7619, 7667, 7782, 7887, 8146, 8386/
data (ix(i),i= 161, 240)/
+ 9048, 9451, 9480, 9825, 10036, 10411, 10488, 10555,
1 8319, 8377, 8540, 8070, 7724, 6976, 6659, 6045,
2 6045, 5911, 5920, 5604, 5546, 5268, 5239, 4961,
3 4846, 4826, 4970, 4855, 4654, 4529, 4270, 4318,
4 4289, 32162, 36452, 36576, 36288, 36010, 35751, 35703,
5 35924, 36144, 36039, 35962, 35761, 35598, 35665, 35598,
6 35780, 35655, 35674, 35540, 35559, 35722, 35636, 35626,
7 35684, 35684, 35463, 35588, 32441, 32565, 32498, 32681,
8 32662, 33055, 33276, 33755, 34034, 34542, 34907, 35741,
9 36116, 42477, 41988, 41614, 41220, 40942, 40549, 39464/
data (ix(i),i= 241, 320)/
+ 39186, 38649, 38227, 37622, 36730, 36336, 36260, 35991,
1 35828, 36269, 36308, 36317, 36423, 36471, 36308, 35790,
2 22107, 21877, 21570, 21378, 21426, 20361, 19843, 19660,
3 19171, 18893, 18586, 18240, 18211, 18077, 18077, 17617,
4 17281, 16954, 16801, 16686, 16081, 15995, 15995, 15822,
5 15938, 16187, 22663, 22318, 22184, 22270, 22596, 22702,
6 22769, 22376, 22021, 22059, 22280, 22347, 22078, 23738,
7 22568, 22232, 21915, 21982, 22193, 23153, 28075, 34600,
8 34139, 33794, 33122, 32345, 31942, 31529, 29965, 28939,
9 28468, 28219, 25197, 25705, 27912, 24026, 24986, 35166/
data (ix(i),i= 321, 400)/
+ 35626, 36010, 37210, 37305, 37651, 38236, 38601, 40232,
1 41028, 41230, 41441, 41412, 41719, 42065, 42151, 42708,
2 43293, 44022, 44166, 44463, 45020, 45336, 48157, 48186,
3 47524, 48724, 48618, 48762, 48378, 48378, 48244, 47467,
4 47217, 47505, 47476, 47064, 46488, 46421, 46190, 46286,
5 46171, 50921, 48954, 49846, 50028, 51708, 52072, 53713,
6 47476, 48129, 48752, 48915, 49290, 49318, 46968, 46152,
7 46430, 46363, 46325, 46325, 45317, 45116, 43801, 43111,
8 54126, 54241, 54289, 54193, 54442, 54701, 54999, 54932,
9 54999, 54845, 54164, 53991, 53991, 53780, 53665, 53627/
data (ix(i),i= 401, 480)/
+ 53780, 54087, 54173, 54116, 53838, 53751, 53790, 54068,
1 54097, 53982, 53608, 53579, 53511, 53521, 53224, 53271,
2 53204, 53157, 53089, 52792, 52408, 52360, 51938, 51554,
3 51381, 51362, 51055, 7993, 7580, 7532, 7667, 17281,
4 17434, 17492, 17665, 17559, 17396, 17434, 17511, 17415,
5 16820, 16609, 16830, 17175, 17290, 17300, 17549, 16964,
6 25168, 25600, 24726, 24112, 23853, 23134, 22798, 22213,
7 22078, 21618, 21455, 21253, 21052, 20936, 20975, 20821,
8 20754, 20457, 19833, 19344, 19065, 18653, 18576, 18365,
9 18346, 18240, 17847, 26233, 26089, 26099, 25868, 25907/
data (ix(i),i= 481, 560)/
+ 26252, 26252, 26952, 27240, 27605, 28344, 28977, 29025,
1 29409, 29716, 30196, 30378, 30445, 30531, 30397, 30435,
2 30992, 31088, 30934, 30426, 30359, 30071, 29917, 29754,
3 30368, 30503, 30752, 30848, 30733, 30551, 30426, 30474,
4 31002, 31126, 31405, 31702, 31894, 31942, 32172, 32201,
5 32211, 32441, 32479, 32786, 32815, 33093, 33228, 33391,
6 33688, 34034, 34149, 34158, 34369, 34197, 33621, 33458,
7 33295, 33208, 32978, 32450, 32335, 32518, 32489, 32690,
8 32383, 32105, 31990, 31990, 31760, 31654, 31424, 31194,
9 31050, 31069, 30829, 30848, 34859, 34676, 34312, 33784/
data (ix(i),i= 561, 640)/
+ 33160, 32930, 32748, 32537, 31740, 31952, 31952, 34629,
1 34619, 34369, 34484, 34283, 34417, 34446, 34609, 34724,
2 34878, 35079, 34964, 34638, 34475, 34542, 34465, 34513,
3 34580, 35060, 35214, 35329, 35291, 35454, 35578, 35722,
4 35838, 35857, 36192, 36125, 35847, 35847, 36125, 36471,
5 36835, 36902, 36663, 36663, 36432, 36691, 37066, 37104,
6 36787, 36931, 36701, 36701, 37075, 37114, 37459, 37555,
7 37737, 37881, 37900, 38198, 38783, 39311, 39455, 39320,
8 39368, 39628, 40501, 40587, 40568, 40894, 40971, 41633,
9 43072, 43379, 43475, 43379, 43379, 43082, 42813, 42583/
data (ix(i),i= 641, 720)/
+ 42612, 34264, 33919, 33755, 33477, 33343, 33132, 32738,
1 32374, 32191, 31107, 30762, 30551, 30081, 29831, 29601,
2 29198, 29457, 29447, 29543, 29447, 29198, 28871, 28785,
3 28018, 27989, 27874, 28075, 28209, 28881, 29629, 29879,
4 27499, 27451, 26147, 25235, 25197, 24803, 24650, 24563,
5 24381, 24362, 24986, 24870, 24995, 25187, 25072, 25168,
6 25398, 25494, 25475, 25590, 25955, 26089, 26415, 26434,
7 25686, 25734, 25648, 25369, 25350, 25513, 26329, 26463,
8 26501, 26895, 27192, 27240, 27077, 26924, 27116, 27269,
9 27317, 27231, 27260, 27538, 27768, 28229, 28286, 28526/
data (ix(i),i= 721, 800)/
+ 28564, 29188, 24199, 24016, 23738, 23738, 23489, 23278,
1 23191, 23249, 23335, 23067, 22443, 21973, 21925, 21685,
2 21618, 21080, 20783, 20572, 20428, 19919, 19689, 22347,
3 22376, 22165, 22337, 22289, 24237, 24371, 24621, 24544,
4 24890, 24544, 24448, 24343, 24323, 23834, 23489, 22875,
5 24026, 24208, 24640, 24832, 25091, 25388, 25638, 25705,
6 25830, 26051, 26147, 26425, 26655, 27231, 27317, 27250,
7 27048, 26876, 27010, 26962, 26665, 26501, 26482, 26607,
8 26933, 27164, 27557, 28200, 28363, 28209, 28622, 28891,
9 29083, 29006, 29121, 29092, 27950, 27595, 27317, 32326/
data (ix(i),i= 801, 880)/
+ 32230, 31884, 31808, 31597, 31069, 30819, 30982, 31309,
1 31155, 30954, 30675, 30397, 30215, 29927, 29543, 29418,
2 29303, 29121, 28910, 28612, 28430, 28411, 31654, 31932,
3 31952, 32162, 32211, 32211, 32565, 32633, 32633, 33007,
4 33084, 33525, 35127, 35223, 35022, 35099, 34820, 34753,
5 34302, 34312, 34197, 33871, 33717, 33717, 33602, 33208,
6 33064, 37632, 37200, 36807, 36480, 36317, 36087, 35972,
7 35905, 35933, 35770, 35540, 34916, 34705, 34178, 33669,
8 33650, 33871, 34264, 34197, 34283, 34235, 34293, 34302,
9 34321, 34571, 34897, 35147, 35406, 35588, 36116, 36029/
data (ix(i),i= 881, 960)/
+ 35703, 35694, 35559, 34273, 33256, 33103, 33218, 33496,
1 33746, 33679, 33554, 32738, 32930, 33113, 33276, 33391,
2 33909, 33055, 33064, 32441, 32162, 31865, 31443, 31299,
3 30973, 30877, 31462, 31625, 31664, 31760, 28190, 28651,
4 28574, 28785, 29313, 29908, 30215, 30301, 30723, 30531,
5 30579, 30627, 30714, 34859, 35070, 35166, 35463, 35358,
6 35079, 35012, 35041, 25513, 42775, 42612, 42717, 42928,
7 43014, 45519, 45864, 45893, 46008, 46085, 46546, 46891,
8 47073, 47467, 48090, 47179, 47035, 46910, 46757, 46728,
9 46450, 46267, 46037, 45605, 44991, 44531, 44569, 44060/
data (ix(i),i= 961, 968)/
+ 43840, 50787, 50700, 49558, 49376, 48896, 53723, 53281/
data (iy(i),i= 1, 80)/
+ 46567, 46452, 45905, 45454, 45080, 44763, 44715, 44437,
1 44188, 43842, 43679, 43957, 44159, 44236, 43871, 43622,
2 43506, 43363, 43017, 42672, 42441, 41799, 41626, 41184,
3 40570, 40263, 39688, 39304, 38997, 39102, 39333, 39285,
4 38757, 38690, 38344, 38009, 37500, 37020, 36752, 36550,
5 36521, 36819, 37769, 38037, 38239, 39198, 39899, 39899,
6 40139, 40062, 40340, 40388, 40849, 42115, 43113, 43228,
7 44370, 44303, 44293, 44735, 44802, 45032, 45272, 45464,
8 46375, 46414, 46289, 46404, 46183, 46135, 46020, 46164,
9 46231, 46279, 46529, 46711, 53850, 55299, 55145, 52583/
data (iy(i),i= 81, 160)/
+ 52420, 42950, 43363, 43737, 43670, 43785, 44341, 44562,
1 44648, 44994, 45608, 46011, 46116, 46068, 46020, 46183,
2 46471, 47009, 47488, 47690, 47920, 48419, 48841, 48995,
3 49225, 50233, 50693, 51134, 51710, 52257, 52487, 52698,
4 52718, 53312, 53428, 53274, 53668, 47690, 52382, 52564,
5 53169, 53111, 52622, 52276, 51662, 51269, 50770, 50223,
6 49945, 49369, 48774, 48640, 48314, 48285, 48006, 47959,
7 45368, 46039, 46279, 46491, 47843, 48525, 49139, 49158,
8 49340, 49148, 49129, 48937, 48381, 47719, 55798, 55692,
9 53879, 53907, 52718, 52439, 52228, 50683, 50194, 49964/
data (iy(i),i= 161, 240)/
+ 49235, 48544, 48333, 48218, 47469, 47239, 46845, 46682,
1 46922, 47268, 47613, 48237, 48640, 49494, 49734, 50290,
2 50472, 50664, 51298, 51614, 52420, 52651, 53111, 53572,
3 53792, 54023, 54310, 54531, 54646, 55107, 55567, 55750,
4 55807, 45205, 42058, 41942, 40935, 40388, 39189, 39208,
5 40388, 40858, 40935, 40781, 40849, 40686, 40522, 40206,
6 40043, 39832, 39467, 39237, 39169, 39074, 38661, 37529,
7 36560, 36147, 35667, 35188, 45099, 44898, 44639, 44293,
8 44063, 43353, 42825, 42720, 42489, 42403, 42153, 42086,
9 41885, 27646, 27406, 27550, 28135, 28202, 28452, 28509/
data (iy(i),i= 241, 320)/
+ 28740, 29383, 30323, 31148, 32386, 33336, 33432, 33729,
1 34276, 33854, 33499, 33422, 33432, 33681, 33892, 34420,
2 35552, 35600, 36061, 36243, 36358, 37481, 38296, 38623,
3 38862, 39419, 39419, 39784, 40043, 40302, 40446, 40743,
4 40993, 41309, 41472, 41309, 41645, 41664, 41827, 42173,
5 42499, 42806, 32060, 32290, 32328, 32510, 32779, 33096,
6 33556, 33863, 34698, 34890, 34689, 34756, 35159, 29200,
7 31829, 31733, 31618, 31407, 30812, 29488, 26869, 23846,
8 24249, 24249, 24422, 24624, 24988, 25036, 25641, 26389,
9 26869, 26926, 28222, 27598, 26811, 29056, 28490, 23654/
data (iy(i),i= 321, 400)/
+ 23367, 23232, 23194, 23031, 22877, 22647, 22705, 23338,
1 24134, 24220, 24134, 24057, 23885, 23865, 24000, 23827,
2 28039, 28193, 28116, 28260, 28471, 28682, 27329, 26312,
3 26216, 33096, 32280, 31349, 30918, 30342, 30045, 29411,
4 29162, 28864, 28768, 28452, 28346, 28481, 28557, 28759,
5 28826, 27540, 33489, 33979, 34007, 34449, 34717, 34650,
6 26149, 25583, 25218, 24873, 24691, 24067, 23885, 22340,
7 21946, 21783, 21275, 20910, 21754, 22023, 23011, 23597,
8 34660, 34756, 34804, 34890, 34986, 34852, 34487, 34286,
9 34075, 33566, 32664, 31925, 31532, 31349, 31340, 31225/
data (iy(i),i= 401, 480)/
+ 31071, 31167, 31119, 30956, 30745, 30553, 30496, 30659,
1 30544, 30169, 28740, 28298, 28308, 28730, 29056, 29498,
2 29613, 29565, 29402, 29104, 28816, 28615, 27809, 27953,
3 27799, 27521, 27492, 53965, 54301, 54531, 55040, 52094,
4 52027, 51058, 50069, 48851, 48669, 47469, 46615, 46519,
5 46721, 46491, 45867, 45301, 45099, 44408, 44101, 43516,
6 45627, 43353, 43526, 42739, 42873, 43075, 42998, 43391,
7 43334, 43679, 43564, 42988, 42864, 42748, 42336, 42221,
8 41875, 41578, 41472, 42441, 42499, 42662, 43334, 43679,
9 43929, 44072, 44245, 42595, 42086, 41530, 41002, 40590/
data (iy(i),i= 481, 560)/
+ 40225, 40033, 39678, 40273, 40426, 40321, 39880, 39745,
1 39755, 39726, 39294, 39323, 39812, 40043, 40177, 40359,
2 40570, 40734, 41002, 41511, 41952, 42307, 42403, 42720,
3 43286, 43219, 43497, 43842, 43996, 43996, 44207, 44389,
4 44687, 45109, 45320, 45090, 44859, 44677, 44581, 44197,
5 43737, 43602, 43324, 43209, 42777, 42633, 42288, 42259,
6 42067, 42163, 42048, 41453, 41060, 40781, 40302, 40369,
7 40330, 39784, 39736, 39438, 39265, 39093, 38536, 38143,
8 37817, 37769, 37605, 37174, 37174, 36828, 36848, 36742,
9 36838, 37740, 38057, 38584, 35399, 35284, 35303, 35188/
data (iy(i),i= 561, 640)/
+ 35427, 35773, 35821, 35686, 36061, 36368, 36502, 34967,
1 34823, 34497, 34180, 33873, 33710, 33393, 33029, 32962,
2 33010, 32846, 32367, 32367, 32204, 31992, 31714, 31561,
3 31484, 31858, 32021, 31906, 31647, 31743, 32261, 32328,
4 32213, 32069, 31868, 31685, 31637, 31311, 31110, 31292,
5 31119, 31043, 30697, 30285, 30045, 29661, 29430, 29363,
6 29200, 28701, 28653, 28423, 28241, 28049, 28097, 28135,
7 27972, 28135, 28270, 28174, 27473, 27358, 27032, 26667,
8 26485, 26235, 26610, 26523, 26264, 26043, 25784, 25708,
9 25717, 26063, 26197, 26341, 26379, 26773, 26917, 27109/
data (iy(i),i= 641, 720)/
+ 27301, 33019, 32971, 33115, 33230, 33566, 33796, 33470,
1 33624, 33460, 33844, 33547, 33633, 34017, 34017, 34199,
2 34238, 34573, 34919, 35399, 35629, 35648, 35370, 35370,
3 36166, 36531, 36761, 37059, 37088, 37605, 38277, 38651,
4 39553, 39160, 39074, 38392, 37605, 37193, 36684, 36550,
5 35255, 35120, 35351, 35610, 35926, 35792, 35168, 35053,
6 35591, 35504, 34737, 34669, 35130, 35360, 35053, 34890,
7 34094, 33892, 33662, 33566, 33058, 33010, 32770, 32770,
8 33134, 33163, 33412, 33624, 33806, 34314, 34497, 34516,
9 34660, 34890, 35025, 35130, 35360, 34986, 34804, 34487/
data (iy(i),i= 721, 800)/
+ 34535, 34238, 35370, 35293, 35571, 35667, 35715, 35667,
1 35706, 35907, 36176, 36435, 36598, 37116, 37538, 38037,
2 38613, 39093, 38910, 39045, 39553, 40052, 40781, 36483,
3 36281, 36061, 35773, 35629, 34698, 34410, 34113, 33652,
4 33307, 32981, 32962, 32194, 32156, 32530, 32722, 32501,
5 29498, 29498, 29853, 29709, 29642, 29814, 29555, 28749,
6 28797, 28980, 28912, 28874, 29200, 29363, 29661, 29757,
7 29776, 30169, 30486, 30582, 30601, 30678, 30812, 30995,
8 30995, 31158, 31302, 31513, 31762, 32252, 32885, 33067,
9 33249, 33547, 33710, 34103, 33835, 33931, 34017, 33163/
data (iy(i),i= 801, 880)/
+ 33000, 32952, 32818, 32539, 32482, 32050, 31839, 31033,
1 30735, 30572, 30688, 30591, 30774, 30908, 30793, 30860,
2 31225, 31359, 31177, 31129, 31196, 31359, 30640, 30995,
3 31139, 31206, 31417, 31685, 31887, 31973, 32232, 32597,
4 32942, 33000, 31369, 31004, 30659, 30217, 30179, 30083,
5 30054, 30189, 30438, 30438, 30486, 30764, 30831, 30735,
6 30419, 27646, 27339, 27358, 27080, 27147, 27608, 27579,
7 27281, 27032, 26965, 27080, 26840, 26993, 27013, 27473,
8 27655, 27867, 27934, 28250, 28394, 28577, 28797, 28855,
9 29728, 29670, 29277, 29095, 29277, 29306, 29536, 29632/
data (iy(i),i= 881, 960)/
+ 29766, 29958, 30093, 29795, 29997, 29392, 29210, 29056,
1 29056, 29613, 29862, 28318, 28010, 27895, 28039, 28010,
2 28922, 29008, 28749, 28385, 28145, 28145, 27770, 27790,
3 28548, 28596, 29527, 29728, 30381, 30515, 27320, 27511,
4 27857, 28001, 27790, 27943, 27780, 27828, 27761, 28366,
5 28596, 28548, 28481, 26072, 25641, 25458, 25247, 24576,
6 24441, 24192, 24077, 28682, 24009, 24393, 24883, 25190,
7 25650, 28471, 28462, 28462, 28116, 27790, 27473, 27502,
8 27732, 27665, 27579, 26581, 26581, 26878, 26993, 27157,
9 27291, 27291, 27080, 26888, 26264, 26763, 27176, 27358/
data (iy(i),i= 961, 968)/
+ 26715, 29268, 31014, 32309, 32357, 32520, 33710, 33000/
c
data (b1(i),i= 1, 100)/
+ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10,
1 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
2 21, 22, 23, 24, 25, 26, 27, 28, 29, 30,
3 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
4 41, 42, 43, 44, 45, 46, 47, 48, 49, 50,
5 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
6 61, 62, 63, 64, 65, 66, 67, 68, 69, 70,
7 71, 72, 73, 74, 75, 78, 79, 80, 82, 83,
8 84, 85, 86, 87, 88, 89, 90, 91, 92, 93,
9 94, 95, 96, 97, 98, 99, 100, 101, 102, 103/
data (b1(i),i= 101, 200)/
+ 104, 105, 106, 107, 108, 109, 110, 111, 112, 113,
1 114, 115, 116, 77, 81, 119, 120, 121, 122, 123,
2 124, 125, 126, 127, 128, 129, 130, 131, 132, 133,
3 134, 135, 118, 137, 138, 139, 140, 141, 142, 143,
4 144, 145, 146, 147, 148, 149, 118, 151, 78, 77,
5 153, 154, 155, 156, 157, 158, 159, 160, 161, 162,
6 163, 164, 165, 166, 167, 1, 76, 169, 170, 171,
7 172, 173, 174, 175, 176, 177, 178, 179, 180, 181,
8 182, 183, 184, 185, 186, 187, 188, 189, 190, 191,
9 192, 151, 137, 195, 196, 197, 198, 199, 200, 201/
data (b1(i),i= 201, 300)/
+ 202, 203, 204, 205, 206, 207, 208, 209, 210, 211,
1 212, 213, 214, 215, 216, 217, 218, 219, 194, 221,
2 222, 223, 224, 225, 226, 227, 228, 229, 230, 231,
3 232, 195, 234, 235, 236, 237, 238, 239, 240, 241,
4 242, 243, 244, 245, 246, 247, 248, 249, 250, 251,
5 252, 253, 254, 255, 220, 257, 258, 259, 260, 261,
6 262, 263, 264, 265, 266, 267, 268, 269, 270, 271,
7 272, 273, 274, 275, 276, 277, 278, 279, 280, 281,
8 82, 283, 284, 285, 286, 287, 288, 289, 290, 291,
9 292, 293, 294, 257, 283, 297, 298, 299, 300, 301/
data (b1(i),i= 301, 400)/
+ 296, 304, 305, 306, 307, 308, 309, 310, 311, 312,
1 313, 303, 315, 316, 303, 296, 318, 315, 304, 320,
2 321, 322, 323, 324, 325, 326, 327, 328, 329, 330,
3 331, 332, 333, 334, 335, 234, 337, 338, 339, 340,
4 341, 343, 344, 346, 347, 348, 349, 350, 351, 352,
5 353, 354, 355, 356, 357, 358, 359, 360, 342, 343,
6 346, 363, 364, 365, 366, 367, 345, 369, 370, 371,
7 372, 373, 374, 375, 376, 377, 378, 379, 380, 381,
8 382, 383, 336, 368, 385, 386, 387, 388, 389, 390,
9 391, 392, 393, 394, 395, 396, 397, 398, 399, 400/
data (b1(i),i= 401, 500)/
+ 401, 402, 403, 404, 405, 406, 407, 408, 409, 410,
1 411, 412, 413, 414, 415, 416, 417, 418, 419, 420,
2 421, 422, 423, 424, 425, 426, 362, 1, 77, 428,
3 429, 430, 431, 81, 432, 433, 434, 435, 436, 437,
4 438, 439, 440, 441, 442, 443, 444, 445, 446, 447,
5 448, 118, 449, 450, 451, 452, 453, 454, 455, 456,
6 457, 458, 459, 460, 461, 462, 463, 464, 465, 466,
7 467, 468, 469, 470, 471, 472, 473, 474, 447, 450,
8 476, 477, 478, 479, 480, 481, 482, 483, 484, 485,
9 486, 487, 488, 489, 490, 491, 492, 493, 494, 495/
data (b1(i),i= 501, 600)/
+ 496, 497, 498, 499, 500, 501, 502, 503, 504, 505,
1 506, 507, 508, 509, 510, 511, 512, 513, 514, 515,
2 137, 194, 517, 518, 519, 520, 521, 522, 523, 524,
3 525, 526, 527, 528, 529, 530, 531, 532, 533, 534,
4 535, 536, 537, 538, 539, 540, 541, 542, 543, 544,
5 545, 546, 547, 548, 549, 550, 551, 552, 553, 554,
6 555, 492, 220, 557, 558, 559, 560, 561, 562, 563,
7 564, 565, 566, 550, 559, 568, 569, 570, 571, 572,
8 573, 574, 575, 576, 577, 578, 579, 580, 581, 582,
9 583, 584, 585, 586, 587, 588, 589, 590, 591, 592/
data (b1(i),i= 601, 700)/
+ 593, 594, 595, 596, 597, 598, 599, 600, 601, 602,
1 603, 604, 605, 606, 607, 608, 609, 610, 611, 612,
2 613, 614, 615, 616, 617, 618, 619, 620, 621, 622,
3 623, 624, 625, 626, 627, 628, 629, 630, 631, 632,
4 633, 634, 635, 636, 637, 638, 639, 640, 234, 574,
5 642, 643, 644, 645, 646, 647, 648, 649, 650, 651,
6 652, 653, 654, 655, 656, 657, 658, 659, 660, 661,
7 662, 663, 664, 665, 666, 667, 668, 669, 670, 671,
8 492, 484, 673, 674, 675, 676, 677, 678, 679, 680,
9 681, 682, 683, 684, 685, 686, 687, 688, 689, 690/
data (b1(i),i= 701, 800)/
+ 691, 692, 693, 694, 695, 696, 697, 698, 699, 700,
1 701, 702, 703, 704, 705, 706, 707, 708, 709, 710,
2 711, 712, 713, 714, 715, 716, 717, 718, 719, 720,
3 721, 657, 681, 723, 724, 725, 726, 727, 728, 729,
4 730, 731, 732, 733, 734, 735, 736, 737, 738, 739,
5 740, 741, 742, 467, 733, 744, 745, 746, 747, 257,
6 682, 749, 750, 751, 752, 753, 754, 755, 756, 757,
7 758, 759, 760, 296, 761, 762, 763, 764, 765, 766,
8 767, 768, 769, 770, 771, 772, 773, 774, 775, 776,
9 777, 778, 779, 780, 781, 782, 783, 784, 785, 786/
data (b1(i),i= 801, 900)/
+ 787, 788, 789, 790, 791, 792, 793, 794, 795, 657,
1 720, 797, 798, 712, 650, 800, 801, 802, 803, 804,
2 805, 806, 807, 808, 809, 810, 811, 812, 813, 814,
3 815, 816, 817, 818, 819, 820, 821, 788, 823, 824,
4 825, 826, 827, 828, 829, 830, 831, 832, 833, 645,
5 808, 589, 835, 836, 837, 838, 839, 840, 841, 842,
6 843, 844, 845, 846, 847, 848, 826, 617, 850, 851,
7 852, 853, 854, 855, 856, 857, 858, 859, 860, 861,
8 862, 863, 864, 865, 866, 867, 868, 869, 870, 871,
9 872, 873, 874, 875, 876, 877, 878, 879, 880, 881/
data (b1(i),i= 901, 1000)/
+ 882, 838, 873, 841, 885, 886, 887, 888, 889, 890,
1 885, 892, 893, 894, 895, 865, 871, 889, 887, 898,
2 892, 892, 900, 901, 902, 903, 904, 905, 906, 907,
3 908, 909, 823, 303, 911, 912, 913, 914, 915, 916,
4 917, 918, 919, 920, 921, 922, 906, 861, 924, 925,
5 926, 927, 928, 929, 930, 931, 768, 932, 336, 933,
6 934, 935, 936, 633, 342, 938, 939, 940, 941, 942,
7 943, 944, 945, 946, 947, 345, 948, 949, 950, 951,
8 952, 953, 954, 955, 956, 957, 958, 959, 960, 635,
9 362, 962, 963, 964, 965, 346, 368, 967, 963, 0/
data (b2(i),i= 1, 100)/
+ 2, 3, 4, 5, 6, 7, 8, 9, 10, 11,
1 12, 13, 14, 15, 16, 17, 18, 19, 20, 21,
2 22, 23, 24, 25, 26, 27, 28, 29, 30, 31,
3 32, 33, 34, 35, 36, 37, 38, 39, 40, 41,
4 42, 43, 44, 45, 46, 47, 48, 49, 50, 51,
5 52, 53, 54, 55, 56, 57, 58, 59, 60, 61,
6 62, 63, 64, 65, 66, 67, 68, 69, 70, 71,
7 72, 73, 74, 75, 76, 79, 80, 81, 83, 84,
8 85, 86, 87, 88, 89, 90, 91, 92, 93, 94,
9 95, 96, 97, 98, 99, 100, 101, 102, 103, 104/
data (b2(i),i= 101, 200)/
+ 105, 106, 107, 108, 109, 110, 111, 112, 113, 114,
1 115, 116, 117, 117, 119, 120, 121, 122, 123, 124,
2 125, 126, 127, 128, 129, 130, 131, 132, 133, 134,
3 135, 136, 136, 138, 139, 140, 141, 142, 143, 144,
4 145, 146, 147, 148, 149, 150, 150, 152, 152, 153,
5 154, 155, 156, 157, 158, 159, 160, 161, 162, 163,
6 164, 165, 166, 167, 168, 168, 169, 170, 171, 172,
7 173, 174, 175, 176, 177, 178, 179, 180, 181, 182,
8 183, 184, 185, 186, 187, 188, 189, 190, 191, 192,
9 193, 193, 194, 196, 197, 198, 199, 200, 201, 202/
data (b2(i),i= 201, 300)/
+ 203, 204, 205, 206, 207, 208, 209, 210, 211, 212,
1 213, 214, 215, 216, 217, 218, 219, 220, 221, 222,
2 223, 224, 225, 226, 227, 228, 229, 230, 231, 232,
3 233, 233, 235, 236, 237, 238, 239, 240, 241, 242,
4 243, 244, 245, 246, 247, 248, 249, 250, 251, 252,
5 253, 254, 255, 256, 256, 258, 259, 260, 261, 262,
6 263, 264, 265, 266, 267, 268, 269, 270, 271, 272,
7 273, 274, 275, 276, 277, 278, 279, 280, 281, 282,
8 282, 284, 285, 286, 287, 288, 289, 290, 291, 292,
9 293, 294, 295, 295, 297, 298, 299, 300, 301, 302/
data (b2(i),i= 301, 400)/
+ 302, 305, 306, 307, 308, 309, 310, 311, 312, 313,
1 314, 314, 316, 317, 317, 318, 319, 319, 320, 321,
2 322, 323, 324, 325, 326, 327, 328, 329, 330, 331,
3 332, 333, 334, 335, 336, 337, 338, 339, 340, 341,
4 342, 344, 345, 347, 348, 349, 350, 351, 352, 353,
5 354, 355, 356, 357, 358, 359, 360, 361, 361, 362,
6 363, 364, 365, 366, 367, 368, 369, 370, 371, 372,
7 373, 374, 375, 376, 377, 378, 379, 380, 381, 382,
8 383, 384, 384, 385, 386, 387, 388, 389, 390, 391,
9 392, 393, 394, 395, 396, 397, 398, 399, 400, 401/
data (b2(i),i= 401, 500)/
+ 402, 403, 404, 405, 406, 407, 408, 409, 410, 411,
1 412, 413, 414, 415, 416, 417, 418, 419, 420, 421,
2 422, 423, 424, 425, 426, 427, 427, 76, 428, 429,
3 430, 431, 78, 432, 433, 434, 435, 436, 437, 438,
4 439, 440, 441, 442, 443, 444, 445, 446, 447, 448,
5 82, 449, 450, 451, 452, 453, 454, 455, 456, 457,
6 458, 459, 460, 461, 462, 463, 464, 465, 466, 467,
7 468, 469, 470, 471, 472, 473, 474, 475, 475, 476,
8 477, 478, 479, 480, 481, 482, 483, 484, 485, 486,
9 487, 488, 489, 490, 491, 492, 493, 494, 495, 496/
data (b2(i),i= 501, 600)/
+ 497, 498, 499, 500, 501, 502, 503, 504, 505, 506,
1 507, 508, 509, 510, 511, 512, 513, 514, 515, 516,
2 516, 517, 518, 519, 520, 521, 522, 523, 524, 525,
3 526, 527, 528, 529, 530, 531, 532, 533, 534, 535,
4 536, 537, 538, 539, 540, 541, 542, 543, 544, 545,
5 546, 547, 548, 549, 550, 551, 552, 553, 554, 555,
6 556, 556, 557, 558, 559, 560, 561, 562, 563, 564,
7 565, 566, 567, 567, 568, 569, 570, 571, 572, 573,
8 574, 575, 576, 577, 578, 579, 580, 581, 582, 583,
9 584, 585, 586, 587, 588, 589, 590, 591, 592, 593/
data (b2(i),i= 601, 700)/
+ 594, 595, 596, 597, 598, 599, 600, 601, 602, 603,
1 604, 605, 606, 607, 608, 609, 610, 611, 612, 613,
2 614, 615, 616, 617, 618, 619, 620, 621, 622, 623,
3 624, 625, 626, 627, 628, 629, 630, 631, 632, 633,
4 634, 635, 636, 637, 638, 639, 640, 641, 641, 642,
5 643, 644, 645, 646, 647, 648, 649, 650, 651, 652,
6 653, 654, 655, 656, 657, 658, 659, 660, 661, 662,
7 663, 664, 665, 666, 667, 668, 669, 670, 671, 672,
8 672, 673, 674, 675, 676, 677, 678, 679, 680, 681,
9 682, 683, 684, 685, 686, 687, 688, 689, 690, 691/
data (b2(i),i= 701, 800)/
+ 692, 693, 694, 695, 696, 697, 698, 699, 700, 701,
1 702, 703, 704, 705, 706, 707, 708, 709, 710, 711,
2 712, 713, 714, 715, 716, 717, 718, 719, 720, 721,
3 722, 722, 723, 724, 725, 726, 727, 728, 729, 730,
4 731, 732, 733, 734, 735, 736, 737, 738, 739, 740,
5 741, 742, 743, 743, 744, 745, 746, 747, 748, 748,
6 749, 750, 751, 752, 753, 754, 755, 756, 757, 758,
7 759, 760, 283, 761, 762, 763, 764, 765, 766, 767,
8 768, 769, 770, 771, 772, 773, 774, 775, 776, 777,
9 778, 779, 780, 781, 782, 783, 784, 785, 786, 787/
data (b2(i),i= 801, 900)/
+ 788, 789, 790, 791, 792, 793, 794, 795, 796, 796,
1 797, 798, 799, 799, 800, 801, 802, 803, 804, 805,
2 806, 807, 808, 809, 810, 811, 812, 813, 814, 815,
3 816, 817, 818, 819, 820, 821, 822, 822, 824, 825,
4 826, 827, 828, 829, 830, 831, 832, 833, 834, 834,
5 823, 835, 836, 837, 838, 839, 840, 841, 842, 843,
6 844, 845, 846, 847, 848, 849, 849, 850, 851, 852,
7 853, 854, 855, 856, 857, 858, 859, 860, 861, 862,
8 863, 864, 865, 866, 867, 868, 869, 870, 871, 872,
9 873, 874, 875, 876, 877, 878, 879, 880, 881, 882/
data (b2(i),i= 901, 1000)/
+ 883, 883, 884, 884, 886, 887, 888, 889, 890, 891,
1 891, 893, 894, 895, 896, 896, 897, 897, 898, 899,
2 899, 900, 901, 902, 903, 904, 905, 906, 907, 908,
3 909, 910, 910, 911, 912, 913, 914, 915, 916, 917,
4 918, 919, 920, 921, 922, 923, 923, 924, 925, 926,
5 927, 928, 929, 930, 931, 304, 932, 315, 933, 934,
6 935, 936, 937, 937, 938, 939, 940, 941, 942, 943,
7 944, 945, 946, 947, 343, 948, 949, 950, 951, 952,
8 953, 954, 955, 956, 957, 958, 959, 960, 961, 961,
9 962, 963, 964, 965, 966, 966, 967, 968, 968, 0/
c
sp(2)='mexico'
sp(1)='mexico'
sp(3)='mexico'
c
c istate=0/1 dont use/use the internal edges/vertices that
c define state borders
c
istate=1
c
do i=1,nvf
vx(i)=real(ix(i),rknd)*1.0e-2_rknd
vy(i)=real(iy(i),rknd)*1.0e-2_rknd
enddo
do i=1,nbf
ibndry(1,i)=b1(i)
ibndry(2,i)=b2(i)
ibndry(3,i)=0
ibndry(4,i)=2
if(i>nbound) ibndry(4,i)=0
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=0
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
rp(13)=grade
if(istate==0) then
ip(2)=nbound
ip(3)=nbound
endif
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd20(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 :: ntf,nvf,nbf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), dimension(200) :: xm,ym
real(kind=rknd), save, dimension(5) :: ratio,diam,xcen,ycen
real(kind=rknd), save, dimension(5) :: angle
real(kind=rknd), save :: hmax,grade
character(len=80), dimension(100) :: sp,su
cy
external sxy
data hmax,grade/0.1e0_rknd,1.75e0_rknd/
data ratio/4.0e0_rknd,1.0e0_rknd,
+ 0.2e0_rknd,2.0e0_rknd,0.7e0_rknd/
data diam/4.0e0_rknd,3.0e0_rknd,
+ 4.0e0_rknd,2.0e0_rknd,4.0e0_rknd/
data angle/0.0e0_rknd,0.0e0_rknd,
+ 45.0e0_rknd,30.0e0_rknd,60.0e0_rknd/
data xcen/0.0e0_rknd,2.0e0_rknd,
+ 5.0e0_rknd,7.0e0_rknd,9.0e0_rknd/
data ycen/2.0e0_rknd,-1.0e0_rknd,
+ 3.0e0_rknd,-2.0e0_rknd,1.5e0_rknd/
c
sp(2)='ellipse'
sp(1)='ellipse'
sp(3)='ellipse'
c
ip(1)=0
ip(2)=0
ip(3)=0
c
nel=5
nn=32
c
do i=1,nel-1
ii=i
do j=i+1,nel
if(xcen(j)<xcen(ii)) ii=j
if(xcen(j)==xcen(ii).and.ycen(j)<ycen(ii)) ii=j
enddo
if(ii/=i) then
t=xcen(i)
xcen(i)=xcen(ii)
xcen(ii)=t
t=ycen(i)
ycen(i)=ycen(ii)
ycen(ii)=t
t=diam(i)
diam(i)=diam(ii)
diam(ii)=t
t=ratio(i)
ratio(i)=ratio(ii)
ratio(ii)=t
t=angle(i)
angle(i)=angle(ii)
angle(ii)=t
endif
enddo
c
xmin=xcen(1)
xmax=xmin
ymin=ycen(1)
ymax=ymin
rmax=0.0e0_rknd
do i=1,nel
r=diam(i)/2.0e0_rknd
xmin=min(xcen(i)-r,xmin)
xmax=max(xcen(i)+r,xmax)
ymin=min(ycen(i)-r,ymin)
ymax=max(ycen(i)+r,ymax)
rmax=max(r,rmax)
enddo
c
nbf=0
nvf=0
ncf=0
c
do i=1,nel
call ellpse(vx(nvf+1),vy(nvf+1),xm(ncf+1),ym(ncf+1),
+ ratio(i),angle(i),diam(i),xcen(i),ycen(i),nn)
do k=1,nn
ibndry(1,nbf+k)=nvf+k
ibndry(2,nbf+k)=nvf+k+1
if(k==nn) ibndry(2,nbf+k)=nvf+1
ibndry(3,nbf+k)=ncf+k
ibndry(4,nbf+k)=1
ibndry(5,nbf+k)=0
ibndry(6,nbf+k)=0
ibndry(7,nbf+k)=-(nbf+k)
enddo
nvf=nvf+nn
nbf=nbf+nn
ncf=ncf+nn
enddo
c
vx(nvf+1)=xmin-rmax
vy(nvf+1)=ymin-rmax
vx(nvf+2)=xmax+rmax
vy(nvf+2)=ymin-rmax
vx(nvf+3)=xmax+rmax
vy(nvf+3)=ymax+rmax
vx(nvf+4)=xmin-rmax
vy(nvf+4)=ymax+rmax
c
do i=1,4
ibndry(1,nbf+i)=nvf+i
ibndry(2,nbf+i)=nvf+i+1
ibndry(3,nbf+i)=0
ibndry(4,nbf+i)=2
ibndry(5,nbf+i)=0
ibndry(6,nbf+i)=0
ibndry(7,nbf+i)=nbf+i
enddo
ibndry(2,nbf+4)=nvf+1
nvf=nvf+4
nbf=nbf+4
c
do i=1,nel+1
do j=1,7
ibndry(j,nbf+i)=0
enddo
enddo
c
ii=nvf-3
i=1
k1=i
k2=k1+nn-1
d0=(vx(ii)-vx(i))**2+(vy(ii)-vy(i))**2
do k=k1,k2
dd=(vx(ii)-vx(k))**2+(vy(ii)-vy(k))**2
if(dd<d0) then
d0=dd
i=k
endif
enddo
ibndry(1,nbf+1)=ii
ibndry(2,nbf+1)=i
nbf=nbf+1
c
do m=1,nel-1
ii=1+(m-1)*nn
i=1+m*nn
d0=(vx(ii)-vx(i))**2+(vy(ii)-vy(i))**2
do mm=1,2
k1=1+m*nn
k2=k1+nn-1
do k=k1,k2
dd=(vx(ii)-vx(k))**2+(vy(ii)-vy(k))**2
if(dd<d0) then
d0=dd
i=k
endif
enddo
k1=1+(m-1)*nn
k2=k1+nn-1
do k=k1,k2
dd=(vx(i)-vx(k))**2+(vy(i)-vy(k))**2
if(dd<d0) then
d0=dd
ii=k
endif
enddo
enddo
ibndry(1,nbf+1)=ii
ibndry(2,nbf+1)=i
nbf=nbf+1
enddo
c
ii=nvf-1
i=1+(nel-1)*nn
k1=i
k2=k1+nn-1
d0=(vx(ii)-vx(i))**2+(vy(ii)-vy(i))**2
do k=k1,k2
dd=(vx(ii)-vx(k))**2+(vy(ii)-vy(k))**2
if(dd<d0) then
d0=dd
i=k
endif
enddo
ibndry(1,nbf+1)=ii
ibndry(2,nbf+1)=i
nbf=nbf+1
c
ip(1)=2
ip(2)=nvf
ip(3)=nbf
c
do i=1,nbf
if(ibndry(3,i)>0) then
sf(1,i)=xm(ibndry(3,i))
sf(2,i)=ym(ibndry(3,i))
else
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
endif
enddo
c
c make itnode, find symmetries
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
c
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine ellpse(vx,vy,xm,ym,ratio,angle,diam,xcen,ycen,nn)
cx
use mthdef
implicit real(kind=rknd) (a-h,o-z)
implicit integer(kind=iknd) (i-n)
real(kind=rknd), dimension(*) :: vx,vy,xm,ym
real(kind=rknd), dimension(1025) :: x,y
cy
if(nn>1024) stop 5567
pi=3.141592653589793e0_rknd
if(ratio<1.0e0_rknd) then
rat=1.0e0_rknd/ratio
ang=(angle+90.0e0_rknd)/180.0e0_rknd
else
rat=ratio
ang=angle/180.0e0_rknd
endif
c
c=cos(ang*pi)
s=sin(ang*pi)
do i=1,2*nn+1
ang=real(i-1,rknd)*pi/real(nn,rknd)
xx=cos(ang)/2.0e0_rknd
yy=sin(ang)/(2.0e0_rknd*rat)
x(i)=(c*xx+s*yy)*diam+xcen
y(i)=(c*yy-s*xx)*diam+ycen
enddo
do i=1,nn
ii=2*i
vx(i)=x(ii-1)
vy(i)=y(ii-1)
call centre(x(ii-1),y(ii-1),x(ii+1),y(ii+1),
+ x(ii),y(ii),xm(i),ym(i))
enddo
c
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine fixit(n,x,y)
cx
use mthdef
implicit real(kind=rknd) (a-h,o-z)
implicit integer(kind=iknd) (i-n)
integer(kind=iknd), dimension(20000) :: ix,iy
real(kind=rknd), dimension(*) :: x,y
character(len=1), dimension(10) :: mark
cy
data mark/'+','1','2','3','4','5','6','7','8','9'/
scale=1.0e1_rknd
do i=1,n
if(x(i)>0.0e0_rknd) then
ix(i)=int(x(i)*scale+0.5e0_rknd)
else if(x(i)<0) then
ix(i)=int(x(i)*scale-0.5e0_rknd)
else
ix(i)=0
endif
if(y(i)>0) then
iy(i)=int(y(i)*scale+0.5e0_rknd)
else if(y(i)<0) then
iy(i)=int(y(i)*scale-0.5e0_rknd)
else
iy(i)=0
endif
enddo
do i=n+1,n+100
ix(i)=0
iy(i)=0
enddo
101 format(12x,'data (ix(i),i=',i4,',',i4,')/')
104 format(12x,'data (iy(i),i=',i4,',',i4,')/')
102 format(5x,a1,5x,12(i4,','))
103 format(5x,a1,5x,11(i4,','),i4,'/')
nlin=12
do i=1,n,10*nlin
write(8,101) i,min(i+10*nlin-1,n)
do j=1,10
istrt=i+(j-1)*nlin
if(j/=10) then
write(8,102) mark(j),(ix(k),k=istrt,istrt+nlin-1)
else
write(8,103) mark(j),(ix(k),k=istrt,istrt+nlin-1)
endif
enddo
enddo
do i=1,n,10*nlin
write(8,104) i,min(i+10*nlin-1,n)
do j=1,10
istrt=i+(j-1)*nlin
if(j/=10) then
write(8,102) mark(j),(iy(k),k=istrt,istrt+nlin-1)
else
write(8,103) mark(j),(iy(k),k=istrt,istrt+nlin-1)
endif
enddo
enddo
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd21(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 :: ntf,nbf,nvf
real(kind=rknd), dimension(*) :: vx,vy
real(kind=rknd), dimension(2,*) :: sf
real(kind=rknd), dimension(100) :: rp,ru
real(kind=rknd), save :: hmax
character(len=80), dimension(100) :: sp,su
cy
external sxy
data ntf,nbf,nvf/2,12,10/
data hmax/0.5e0_rknd/
c
sp(2)='hole'
sp(1)='hole'
sp(3)='hole'
c
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
rp(12)=hmax
c
vx(1)=-2.0e0_rknd
vy(1)=2.0e0_rknd
vx(2)=-2.0e0_rknd
vy(2)=-2.0e0_rknd
vx(3)=0.0e0_rknd
vy(3)=-2.0e0_rknd
vx(4)=2.0e0_rknd
vy(4)=-2.0e0_rknd
vx(5)=2.0e0_rknd
vy(5)=2.0e0_rknd
vx(6)=0.0e0_rknd
vy(6)=2.0e0_rknd
vx(7)=0.0e0_rknd
vy(7)=1.0e0_rknd
vx(8)=-1.0e0_rknd
vy(8)=0.0e0_rknd
vx(9)=0.0e0_rknd
vy(9)=-1.0e0_rknd
vx(10)=1.0e0_rknd
vy(10)=0.0e0_rknd
c
do i=1,nbf
ibndry(1,i)=i-1
ibndry(2,i)=i
ibndry(3,i)=0
if(i>=8.and.i<=11) ibndry(3,i)=1
ibndry(4,i)=1
if(i>=2.and.i<=5) ibndry(4,i)=2
ibndry(5,i)=0
ibndry(6,i)=0
ibndry(7,i)=0
sf(1,i)=0.0e0_rknd
sf(2,i)=0.0e0_rknd
enddo
c
ibndry(1,1)=6
ibndry(1,12)=3
ibndry(2,11)=7
ibndry(2,12)=9
ibndry(4,7)=0
ibndry(4,12)=0
c
itnode(1,1)=6
itnode(2,1)=1
itnode(3,1)=0
itnode(4,1)=0
itnode(5,1)=1
itnode(1,2)=6
itnode(2,2)=7
itnode(3,2)=-1
itnode(4,2)=0
itnode(5,2)=2
c
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd22(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(2000) :: ix,iy
integer(kind=iknd), save :: ispd=1,iprob=1,itask=0,ifirst=1
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
c
data (ix(i),i= 1, 100)/
+ 55110,55856,56437,56519,56418,56905,57294,57852,58624,59675,
1 60375,60644,60972,61319,61603,61757,61792,61986,62091,62074,
2 61965,61997,62148,62286,62458,62587,62698,62806,62882,63203,
3 63123,63000,62778,62628,62325,62092,61833,61745,61287,61222,
4 61075,60832,60598,59610,58479,58285,57914,57466,56332,55617,
5 55536,55298,55281,55771,55955,56214,56792,56935,56907,56797,
6 56656,56612,56612,56716,56721,56511,56335,56231,56030,55670,
7 55281,55152,55135,54977,54864,54762,54549,54379,54200,53893,
8 53555,53363,53165,52514,52369,52385,52514,52981,53136,53364,
9 53380,53293,53135,53006,52931,52951,52745,52765,52680,52534/
data (ix(i),i= 101, 200)/
+ 52067,51992,51984,52373,52339,52197,52101,51893,51795,51817,
1 51769,51405,51277,51025,50881,50812,50683,50556,50313,50237,
2 50133,50074,49982,49842,49659,49664,49801,49934,50284,50811,
3 51392,51589,51580,51413,50687,50430,50097,50009,49935,49976,
4 50279,50195,50079,49975,49839,49709,49594,49364,49370,49381,
5 49421,49439,49440,49381,49324,49272,49130,49047,48743,48589,
6 48504,48628,48595,48499,48386,48311,48315,48355,48355,48283,
7 48225,48173,48140,48145,48116,48060,47968,47864,47514,47206,
8 47080,47095,47214,47191,47405,47420,47304,47120,46977,46887,
9 46731,46678,46742,46688,46596,46387,46232,46090,45975,46029/
data (ix(i),i= 201, 300)/
+ 45965,45638,45534,45417,44674,44336,44350,44532,44703,44742,
1 44573,44430,44273,44116,44038,43972,44129,44286,44627,45309,
2 45531,45739,45909,45949,45690,45691,45743,45985,46094,46111,
3 46076,45991,45892,45822,45730,45581,45417,45217,45050,44841,
4 44697,44642,44679,44795,44946,45078,45209,45327,45373,45336,
5 45266,45192,44903,44483,44169,44037,43906,43945,43866,43696,
6 43618,43434,42880,42735,42562,42389,42214,42278,42410,43004,
7 43689,44033,46013,46617,47061,47858,49646,49999,50637,51485,
8 52141,52546,54087,54289,54796,28114,28346,28587,28903,29184,
9 29356,29315,29149,28986,28656,28155,28046,28228,28345,28564/
data (ix(i),i= 301, 400)/
+ 28845,28923,28744,28928,29349,29759,29958,30035,30720,30891,
1 30997,31028,31241,31325,31454,31441,31584,31726,31831,31892,
2 31882,31826,31730,31642,31541,31471,31372,31333,31320,31269,
3 31186,31102,30913,30836,30915,30881,30816,30733,29899,29697,
4 29560,29624,29257,28980,28688,28173,27708,27630,27599,27620,
5 27671,27822,27978,28070,28009,27905,27758,27629,27520,27441,
6 27290,27208,27065,26816,26620,26614,26683,26818,26950,27105,
7 27265,27276,27157,26933,26727,26601,26516,26377,26172,26054,
8 25991,25901,25763,25722,25883,26269,26510,26902,27777,27950,
9 27986,20188,24195,28203,32210,36217,40225,44232,48240,52247/
data (ix(i),i= 401, 500)/
+ 56254,60262,60937,61075,64269,64347,64913,65248,65411,65920,
1 65962,66228,66861,67349,68276,72284,76291,76291,76291,76291,
2 76291,76291,76291,76291,76291,76291,76291,76291,76291,76291,
3 76291,76291,76291,76291,76291,76291,76291,76291,76291,76291,
4 76387,76407,76460,76566,76715,76728,76717,76653,76678,76892,
5 76926,76983,77042,77073,77153,77334,77339,77394,77415,77411,
6 77278,77199,77090,76884,76721,76660,76617,76686,76803,76827,
7 76704,76701,76827,76899,77022,77130,77175,77165,77282,77408,
8 77474,77591,77633,77708,77834,77868,77847,77863,77906,78037,
9 78247,78331,78526,78670,78729,78837,78881,78859,78886,78981/
data (ix(i),i= 501, 600)/
+ 79217,79281,79350,79361,79380,79564,79625,79744,79787,79967,
1 80085,80076,80132,80226,80237,80206,80282,80283,80288,80370,
2 80513,80607,80761,80950,81018,81149,81153,81195,81282,81331,
3 81349,81247,81232,81151,81177,81238,81331,81396,81459,81544,
4 81627,81661,81699,81681,81695,81633,81619,81702,81792,81827,
5 81883,81862,81894,81985,82032,82083,82190,82293,82328,82478,
6 82452,82518,82665,82722,82715,82641,82616,82634,82804,82943,
7 82951,82935,83075,83099,83104,82901,82873,82882,82956,82943,
8 82985,83138,83177,83189,83118,83122,83231,83296,83347,83495,
9 83496,83372,83309,83310,83346,83398,83562,83605,83705,83761/
data (ix(i),i= 601, 700)/
+ 83890,83918,84039,84070,84089,84074,84170,84236,84295,84373,
1 84385,84538,84695,84733,84893,85149,85188,85227,85327,85336,
2 85287,85247,85083,84990,84931,84909,84953,85054,85070,85415,
3 85492,85584,85710,85833,85877,85915,86089,86387,86603,86683,
4 86822,86832,86820,86868,86890,87018,87006,86998,87020,87217,
5 87280,87424,87444,87400,87427,87600,87636,87774,87801,87948,
6 88155,88162,88294,88330,88413,88431,88455,88444,88581,88606,
7 88751,88802,88938,89019,89100,89069,89072,89282,89350,89445,
8 89583,89631,89674,89755,89735,89738,89820,89843,89978,89990,
9 89979,90016,90169,90342,90453,90483,90530,90549,90693,90744/
data (ix(i),i= 701, 800)/
+ 90768,90756,90791,90844,90916,91089,91135,91176,91179,91142,
1 91128,91209,91189,91251,91301,91425,91481,91637,91714,91700,
2 91730,91827,91856,91982,92053,92171,92220,92314,92314,92226,
3 92224,92313,92328,92401,92686,92722,92755,93091,93270,93281,
4 93414,93556,93689,93717,93805,93874,93953,94015,94044,93944,
5 93637,93623,93652,93780,93984,94038,94093,94093,94227,94355,
6 94382,94494,94543,94667,94681,94784,94832,94939,95050,95101,
7 95135,95125,95191,95182,95079,95110,95231,95394,95368,95393,
8 95481,95528,95607,95615,95722,95880,96021,96134,96243,96320,
9 96448,96501,96661,96806,96832,96864,96885,96936,97048,97079/
data (ix(i),i= 801, 900)/
+ 97210,97250,97243,97143,97133,97194,97245,97271,97244,97389,
1 97374,97365,97390,97492,97641,97620,97672,97511,97534,97531,
2 97605,97686,97703,97697,97664,97717,97640,97594,97601,97599,
3 97411,97315,97307,97393,97528,97597,97735,97927,97958,97979,
4 97882,97855,97863,97871,97934,97994,98114,98126,98158,98242,
5 98281,98310,98341,98458,98468,98595,98689,98694,98656,98673,
6 98731,98936,99098,99148,99298,99492,99582,99609,99605,99750,
7 99894,99974,100000,99995,99950,99890,98675,97417,96319,94463,
8 92086,88190,86574,83508,80930,79027,78616,77777,77236,76681,
9 74523,72862,71897,71823,70209,68481,67254,67180,66895,65288/
data (ix(i),i= 901, 1000)/
+ 65208,64904,64831,64628,64438,64024,63994,64044,64008,63915,
1 63791,63602,63702,63305,63215,63256,63361,63536,64052,64499,
2 64670,64816,64709,64603,64493,64274,63684,63187,63075,63024,
3 63084,63000,63179,63224,63215,63247,63310,63460,63500,63489,
4 63415,63361,63246,63212,63167,63094,62953,62748,62566,62288,
5 62186,62124,61906,61581,61236,60470,60372,60150,59919,59839,
6 59879,60025,60131,60244,60668,60872,60984,61033,61263,61579,
7 61636,61538,61082,60841,60473,60563,60577,60325,60344,60430,
8 60817,60995,61017,60902,60482,60367,60246,60292,60731,60812,
9 60728,60167,60025,60029,60274,60305,60256,60101,59952,59876/
data (ix(i),i= 1001, 1100)/
+ 59808,59840,59730,59611,59301,58669,58395,57840,57228,56740,
1 56643,56648,56792,57057,57225,57254,57083,57039,57167,57482,
2 57561,57585,57303,57376,57474,57767,58027,58321,58635,58683,
3 58513,58119,57879,57553,57248,56944,56756,56471,56277,56149,
4 56002,55929,55835,55724,55828,56383,56685,56652,56521,56597,
5 56751,56758,56585,56475,56231,56374,56250,55646,55424,55265,
6 55083,54986,54727,54632,54554,54394,54345,54249,54116,53961,
7 53888,53922,54008,54268,54209,54131,53983,53774,53603,53387,
8 53228,53022,52948,53033,52930,52188,51783,51519,51284,51167,
9 51105,51386,52050,52092,51883,51327,51183,51134,51928,52757/
data (ix(i),i= 1101, 1200)/
+ 53277,53629,53694,53644,53701,53944,54030,53967,53780,53701,
1 53678,53811,53725,53411,53354,53382,53124,52804,52646,51406,
2 51355,51533,51933,51987,51871,51476,51310,51071,50690,50532,
3 50092,49947,49871,49889,50155,51027,51212,51267,51229,51032,
4 50556,50320,50033,49846,49761,49362,49203,48966,48838,48759,
5 48413,48196,48116,47917,47469,47153,46441,46229,45979,45966,
6 46073,46550,47279,48407,48685,49120,49266,49334,49181,49372,
7 49333,49160,48930,48676,46883,46352,46007,45941,45689,45609,
8 45622,45569,45263,44958,44691,44823,46116,46997,47465,47492,
9 47372,46677,46197,46077,46023,46278,46238,46051,45636,45076/
data (ix(i),i= 1201, 1300)/
+ 44902,44862,44928,45222,45302,45931,45944,46078,46562,46858,
1 46954,46912,46764,46468,46051,45716,45461,45019,44898,44804,
2 44817,45044,45246,45286,45218,44949,44586,44423,44433,44658,
3 44792,45265,45399,45492,45762,45965,46047,46102,46238,46359,
4 46318,46399,46943,47311,47530,48129,48520,49018,49311,49406,
5 49476,49450,48968,48549,48321,48376,48567,48963,48950,48787,
6 48459,47980,47748,47489,47312,47271,47325,47531,48039,48314,
7 48356,48268,48068,47944,47806,48081,47888,47558,46928,46751,
8 46573,46176,46039,45971,46150,46042,45280,45072,44936,44777,
9 44655,44501,44662,44536,44562,44699,45151,45178,45053,44655/
data (ix(i),i= 1301, 1400)/
+ 44273,44194,44202,44228,44215,44145,44031,43944,43915,43867,
1 43895,43874,43822,43915,43845,43721,43476,43355,43174,42960,
2 42757,42578,42628,42977,43107,43061,43564,43769,44322,44446,
3 44403,43711,43611,43663,43965,44019,43962,43768,43575,43163,
4 43045,42936,42702,42271,42216,42142,42131,41991,41702,41501,
5 41245,40963,40734,40453,40419,40467,40393,40485,40411,40143,
6 39981,39947,40114,40397,40635,40796,41136,41469,41808,42284,
7 42672,43164,43349,43494,43455,43527,43862,44064,44232,44270,
8 44183,43325,43156,42933,42668,42371,42177,41869,41572,41227,
9 40912,40615,40731,40727,40608,40634,41024,41174,41433,41775/
data (ix(i),i= 1401, 1500)/
+ 41890,41891,41793,41568,41451,41320,41193,41089,41083,41256,
1 41168,41028,40351,40139,39547,39355,38978,38646,38623,38773,
2 38856,38684,38459,38322,38165,37831,37577,37313,36951,36328,
3 35525,35469,35650,35510,35035,34802,34444,34078,33923,34003,
4 34302,34175,34187,34039,34063,34175,34476,34743,34986,35562,
5 35803,35974,36020,35942,35677,35050,34604,34596,34717,35224,
6 35548,35754,35759,35582,35534,35664,35947,36201,36419,36560,
7 36445,36151,36070,36104,36512,36980,37345,37615,37778,37967,
8 38005,37732,37294,37013,36831,36795,36730,36554,36553,36691,
9 36794,36784,36150,35646,35500,35460,35595,35750,36077,36127/
data (ix(i),i= 1501, 1600)/
+ 36131,36140,36150,36160,36168,36181,36209,36226,36242,36256,
1 36259,36255,36247,36234,36223,36213,36204,36198,36156,36058,
2 36057,36055,36053,36052,36046,36042,36039,36035,36034,36031,
3 36028,36024,36016,36003,35989,35991,35995,35999,35998,35994,
4 35998,36001,36014,36023,36026,36017,36009,36008,36008,36004,
5 36009,36008,36003,35996,35994,35983,35966,35964,35955,35951,
6 35951,35949,35941,35926,35915,35909,35897,35883,35872,35842,
7 35831,35823,35820,35813,35813,35807,35798,35789,35788,35790,
8 35788,35781,35780,35783,35780,35777,35772,35764,35760,35769,
9 35777,35776,35765,35760,35757,35756,35748,35743,35734,35723/
data (ix(i),i= 1601, 1700)/
+ 35716,35715,35709,35705,35704,35696,35695,35689,35690,35688,
1 35673,35664,35660,35659,35665,35669,35671,35672,35673,35675,
2 35674,35671,35665,35654,35653,35650,35650,35646,35647,35641,
3 35636,35635,35638,35639,35627,35616,35613,35612,35616,35624,
4 35633,35638,35646,35645,35645,35655,35672,35675,35676,35675,
5 35675,35671,35663,35661,35671,35686,35688,35695,35706,35713,
6 35717,35720,35712,35711,35712,35716,35723,35727,35721,35718,
7 35715,35716,35721,35717,35712,35707,35689,35672,35659,35643,
8 35620,35613,35586,35569,35534,35525,35522,35534,35541,35554,
9 35567,35578,35582,35583,35584,35583,35576,35570,35559,35555/
data (ix(i),i= 1701, 1800)/
+ 35537,35528,35525,35527,35524,35523,35519,35521,35528,35531,
1 35533,35534,35532,35553,35765,35776,35781,35786,35788,35788,
2 35812,35828,35836,35842,35851,35854,35867,35872,35892,35908,
3 35916,35927,35933,35937,35947,35945,35949,35960,35968,35983,
4 35993,35991,35992,35993,35998,36012,36026,36041,36048,36059,
5 36079,36092,36127,36164,36132,36113,36121,36124,36136,36156,
6 36200,36150,36091,36073,35972,35806,35233,34840,34772,34512,
7 34414,34346,34049,33856,33723,33635,33249,33234,33083,33033,
8 32969,32885,32685,32145,31940,31917,31877,31405,30342,30003,
9 29871,28847,28743,28560,28689,27694,27703,27698,27918,27919/
data (ix(i),i= 1801, 1900)/
+ 27980,26695,26723,26715,27170,27174,26747,26692,25945,25583,
1 25127,24711,24136,23889,23873,23841,23488,23382,23099,22719,
2 22357,22377,22350,22346,22655,21403,21023,20769,19055,18817,
3 18885,18888,18535,18642,18384,18272,17949,17451,17382,16354,
4 16083,16093,16063,16167,16179,15792,15820,15788,15888,15895,
5 15983,16075,16075,16070,15721,15638,15603,15567,15528,15512,
6 15322,15254,15241,15211,15204,15162,14267,13285,12757,12738,
7 12390,12152,11410,11144,10768,11222,11225,11234,11234,11071,
8 10963,10738,10707,10272,10302,10193,10223,10207,10211, 9851,
9 8910, 8413, 8047, 8045, 7803, 7805, 7386, 7071, 6373, 6088/
data (ix(i),i= 1901, 1920)/
+ 6018, 6021, 6024, 6023, 6135, 6137, 5756, 2686, 1643, 1330,
1 63, 0, 151, 4158, 7251, 8166,12173,16180,16515,16697/
data (iy(i),i= 1, 100)/
+ 14254,12799,12005,11759,11024, 9838, 9605, 9429, 8873, 8310,
1 7667, 6816, 6460, 5152, 4729, 4231, 3624, 3191, 2750, 2405,
2 1966, 1860, 1863, 2109, 3044, 3118, 3055, 2699, 1949, 1087,
3 750, 554, 386, 383, 585, 451, 9, 0, 110, 324,
4 406, 306, 379, 1014, 1604, 2119, 1940, 2093, 3136, 3630,
5 3876, 3942, 4184, 4482, 4677, 5212, 5812, 6613, 7270, 7794,
6 7751, 7599, 7220, 6903, 6696, 6106, 5790, 5728, 5775, 5642,
7 5787, 5708, 5571, 5391, 5449, 5851, 5932, 5787, 5297, 5193,
8 5366, 5754, 5971, 6329, 6663, 6800, 6880, 6703, 6799, 7235,
9 7372, 7445, 7263, 7184, 7222, 7531, 7953, 8262, 8421, 8204/
data (iy(i),i= 101, 200)/
+ 7915, 7953, 8160, 8951, 9142, 9097, 8826, 8697, 8890, 9285,
1 9425, 9199, 9204, 9386, 9254, 8998, 8917, 9008, 9704, 9741,
2 9676, 9214, 9114, 9153, 9675,10018,10530,10783,10961,10907,
3 11057,11221,11428,11452,11170,11179,11311,11382,11591,11745,
4 12147,12477,12532,12467,12041,11959,12015,12694,12841,13060,
5 13228,13349,13443,13485,13444,13370,13296,13354,14033,14105,
6 14013,13667,13537,13421,13355,13405,13509,13660,13831,13933,
7 14009,13953,13841,13614,13452,13259,13072,13006,12996,13226,
8 13572,13795,14085,14411,15146,15283,15423,15220,15257,15413,
9 15399,15279,15106,14901,14714,14580,14651,14859,15169,15546/
data (iy(i),i= 201, 300)/
+ 15633,15380,15398,15537,15097,15219,15442,15647,16058,16384,
1 16574,17038,17107,17005,17040,17212,17743,18017,18238,18250,
2 18093,17713,17694,17848,18538,18709,18829,18980,19051,19178,
3 19310,19362,19375,19334,19294,19322,19392,19455,19516,19522,
4 19458,19374,19272,19196,19242,19248,19159,19077,19010,18878,
5 18807,18782,18631,18531,18326,18327,18584,18739,18859,18705,
6 18139,18191,19117,19236,19768,20041,20744,21001,21088,21175,
7 21536,21587,20629,20157,19669,19105,18304,18311,18067,17951,
8 17693,17409,16856,16089,14616,41722,41841,41790,41127,40838,
9 40473,40179,40151,40345,40511,40202,39903,39368,38743,38279/
data (iy(i),i= 301, 400)/
+ 38213,38047,37744,37123,36142,35946,35514,35039,33849,33705,
1 33780,34021,34170,34141,34011,33651,33381,33360,33338,33235,
2 33130,33063,33047,33072,33097,33119,33025,32931,32847,32784,
3 32752,32844,33193,33051,32799,32643,32554,32583,34107,34316,
4 34615,35013,35195,35792,36114,36455,37140,37306,37563,37667,
5 37779,37864,37882,37973,38112,38183,38202,38123,38035,37944,
6 37851,37851,37975,38535,38974,39260,39352,39449,39484,39535,
7 39672,39785,39886,39919,40008,40097,40093,40046,40037,40079,
8 40150,40294,40332,40459,40575,40675,40848,40861,41218,41385,
9 41459,99737,99737,99737,99737,99737,99737,99737,99737,99737/
data (iy(i),i= 401, 500)/
+ 99737,99737,99742,99742,99737,99738,99741,99743,99744,99745,
1 99745,99745,99745,99743,99737,99737,99737,98948,96158,95674,
2 91982,91179,88901,85941,82620,82078,76769,74062,73553,66780,
3 65504,64576,62635,59930,56946,54731,51323,49193,48387,46734,
4 46723,46580,46534,46537,46465,46220,46159,45970,45881,45926,
5 45812,45821,45868,45865,45698,45480,45385,45213,45164,45031,
6 44906,44926,44953,45080,45138,45146,44978,44829,44658,44569,
7 44325,44286,44161,44215,44130,44007,43868,43744,43416,43134,
8 43118,43040,42949,42971,42924,42966,43102,43139,43142,43015,
9 42949,42931,42906,42942,42910,42732,42592,42470,42427,42164/
data (iy(i),i= 501, 600)/
+ 41670,41332,41354,41493,41506,41435,41434,41324,41247,41215,
1 41010,40839,40847,40944,41146,41213,41399,41697,41821,41849,
2 41806,41690,41559,41336,41279,41103,41001,40932,40872,40810,
3 40737,40469,40440,40209,40236,40157,40253,40291,40250,40105,
4 40147,40127,40003,39935,39823,39668,39584,39415,39331,39114,
5 38996,38897,38845,38839,38887,38989,38958,39028,39007,38789,
6 38699,38563,38572,38540,38462,38176,38101,37902,37711,37572,
7 37532,37424,37168,37101,36896,36724,36618,36515,36402,36263,
8 36193,36154,36124,36052,35875,35780,35591,35455,35416,35392,
9 35337,35169,35014,34959,34765,34733,34700,34708,34810,34833/
data (iy(i),i= 601, 700)/
+ 34804,34784,34858,34987,35188,35292,35410,35540,35586,35518,
1 35391,35177,34932,34886,34805,34723,34748,34718,34458,34425,
2 34174,34126,34044,33950,33834,33610,33391,33131,33113,33039,
3 33089,33058,32762,32490,32459,32476,32665,32606,32615,32632,
4 32772,32810,32985,33039,33019,32463,32183,32098,31906,31682,
5 31600,31551,31516,31375,31291,31082,31091,30956,30860,30603,
6 30426,30212,29932,29886,29659,29554,29330,29230,29016,28917,
7 28788,28762,28760,28745,28894,29041,29127,29424,29615,29558,
8 29633,29640,29600,29419,29290,29211,29093,28979,28811,28761,
9 28662,28458,28240,28049,27831,27621,27572,27584,27493,27427/
data (iy(i),i= 701, 800)/
+ 27219,27111,26955,26850,26765,26675,26619,26406,26335,26107,
1 25985,25717,25643,25473,25392,25296,25277,25261,25213,25044,
2 24998,24923,24925,24844,24806,24663,24629,24382,24380,24246,
3 24170,23916,23750,23751,23396,23403,23451,23399,23102,23056,
4 22888,22859,22711,22684,22547,22297,21989,21921,21828,21724,
5 21580,21512,21420,21321,21331,21296,21126,21032,20893,20809,
6 20740,20613,20587,20465,20344,20242,20255,20372,20284,20186,
7 20037,19961,19789,19720,19587,19509,19277,19146,18987,18863,
8 18686,18691,18869,18930,19047,19065,19134,19068,19106,18938,
9 18539,18370,18041,17873,17804,17561,17501,17497,17423,17266/
data (iy(i),i= 801, 900)/
+ 17211,17170,17109,17008,16924,16736,16685,16482,16363,15755,
1 15641,15393,15104,14828,14611,14452,14015,13561,13430,13235,
2 13013,12765,12542,12190,12112,11817,11617,11385,11313,11078,
3 10929,10796,10602,10432,10297,10265,10246,10339,10316,10256,
4 9724, 9472, 9181, 9116, 8905, 8804, 8781, 8700, 8645, 8600,
5 8551, 8301, 8239, 8060, 7917, 7822, 7728, 7679, 7556, 7387,
6 7318, 7304, 7130, 7125, 7164, 7020, 6903, 6810, 6693, 6232,
7 6061, 5947, 5846, 5730, 5687, 5553, 5557, 5562, 5561, 5560,
8 5537, 5552, 5594, 5588, 5580, 5583, 5584, 5583, 5577, 5581,
9 5584, 5580, 5565, 5564, 5558, 5564, 5537, 5535, 5583, 5591/
data (iy(i),i= 901, 1000)/
+ 5592, 5601, 5679, 6150, 6169, 6021, 5920, 5655, 5614, 5610,
1 5681, 6201, 6675, 6869, 7362, 7929, 7988, 7919, 8041, 8290,
2 8636, 9469, 9618, 9560, 8914, 8624, 8337, 8350, 8412, 8677,
3 9380, 9750,10413,10483,10671,10849,11024,11263,11415,11531,
4 11553,11505,11406,11287,11102,10882,10719,10542,10417,10144,
5 9877, 9381, 9090, 9240, 9546, 9903,10015,10516,10845,11093,
6 11452,11872,11933,11872,10940,10302,10240,10772,11029,11209,
7 11325,11437,11303,11375,12423,12727,13071,13851,13987,14205,
8 14467,14694,14917,15272,15789,15850,15738,15597,15217,14968,
9 14837,14433,14099,13891,13233,13041,12803,12798,13259,13299/
data (iy(i),i= 1001, 1100)/
+ 13218,13026,12879,12474,12463,12096,12357,12706,13695,14209,
1 14405,14577,14620,14396,14368,14469,14791,15018,15389,15868,
2 16294,16603,17155,17409,17593,17763,17746,17915,18308,18546,
3 18575,18136,17997,18018,17883,17368,17173,17259,17736,17830,
4 17700,17446,17348,17579,18020,19108,20003,20660,21133,21473,
5 21774,22032,21974,21825,21133,20160,19495,18548,18078,17983,
6 18045,18325,18425,18792,18830,18735,18324,18140,18060,18136,
7 18347,19257,19563,20015,20276,20315,20184,19678,18515,18215,
8 18119,18249,18459,18765,18873,18422,18473,18397,18422,18478,
9 18652,19038,19443,19596,19638,19384,19424,19649,19963,19999/
data (iy(i),i= 1101, 1200)/
+ 20354,20820,21368,21508,21712,21942,22248,22974,23499,23537,
1 23056,22499,22193,21708,21504,20969,20602,20409,20399,20155,
2 20294,20701,20944,21062,21204,21219,21586,21439,21590,21578,
3 21265,21303,21512,21907,22071,22008,22037,22155,22259,22352,
4 22402,22513,22848,22733,22254,22007,21994,22103,22536,22572,
5 22357,21709,21573,21491,21637,21866,22186,22394,22964,23170,
6 23237,23111,23100,23250,23210,22924,22972,23142,23731,24276,
7 24380,24315,23874,23674,23450,23543,23839,24011,24151,24101,
8 23809,23689,23709,23985,24502,24758,25023,25320,25656,25759,
9 25898,25548,25502,25555,25693,25982,26086,26139,25920,25478/
data (iy(i),i= 1201, 1300)/
+ 25496,25599,25770,25957,26266,26707,27531,27959,28263,28362,
1 28575,28653,28689,28504,28371,28117,27569,27332,27384,27625,
2 27762,27950,27948,28017,28189,28448,28432,28844,29497,30543,
3 30714,30970,31140,31758,32443,32613,32578,32183,32181,32472,
4 32747,32970,33393,33835,33969,33736,32733,32052,31289,31304,
5 31560,31802,32792,33351,34111,34316,34433,34459,34666,34738,
6 34710,34460,34447,34588,35123,35570,35776,36116,36641,37152,
7 37564,38433,38720,38773,38690,37912,37108,36563,35918,35250,
8 35183,35410,35325,35154,34706,34381,33905,34421,34335,33546,
9 33512,34062,34594,34989,35091,35177,35073,35176,35399,35538/
data (iy(i),i= 1301, 1400)/
+ 35298,34989,34795,34705,34503,34307,34227,34301,34416,34556,
1 34654,34740,35228,35503,35675,35726,35571,35364,34385,33990,
2 33817,33919,34296,34865,35466,35912,37613,37785,37615,37649,
3 37838,38094,38334,38626,38816,38936,39073,39124,39089,38829,
4 38365,38210,38191,38394,38528,38701,39387,40295,40978,41284,
5 41555,41722,42010,42742,42982,43360,43530,43806,43976,44022,
6 44311,44550,44656,44576,44598,44309,44177,43651,43518,43558,
7 43802,43720,43533,43190,42950,42779,42814,42335,42266,42593,
8 42884,44168,44235,44096,44059,44177,44054,44034,44150,44539,
9 44774,45455,45920,46092,46364,46467,46715,46940,47039,46975/
data (iy(i),i= 1401, 1500)/
+ 47012,47126,47196,47206,47209,47245,47304,47539,47796,48297,
1 48415,48327,47028,46870,46428,46217,45539,45325,45135,44797,
2 44370,43474,43502,43671,43925,44618,44813,45114,45839,46656,
3 48013,49074,49741,49847,49929,50348,51343,51497,51747,52196,
4 52570,52839,53199,53586,53690,53764,53521,53036,52756,51374,
5 51093,51117,51410,51578,52065,52930,53870,54041,54252,54583,
6 55454,55720,55943,56004,56105,56453,56636,56492,56501,56789,
7 56913,56987,57156,57397,58099,59146,60257,60472,60443,60518,
8 60673,60853,60752,60794,61250,61712,61847,61823,61514,61073,
9 60391,59945,57999,57190,57184,57422,57994,58446,59008,59218/
data (iy(i),i= 1501, 1600)/
+ 59222,59225,59235,59248,59255,59266,59286,59307,59341,59371,
1 59387,59428,59442,59472,59492,59512,59531,59541,59852,60055,
2 60072,60101,60133,60158,60186,60196,60208,60235,60272,60300,
3 60325,60351,60369,60368,60367,60374,60403,60415,60440,60464,
4 60483,60489,60496,60498,60505,60535,60557,60590,60623,60659,
5 60705,60723,60741,60753,60771,60795,60811,60825,60847,60867,
6 60886,60896,60910,60933,60953,60972,60987,61003,61017,61056,
7 61072,61084,61090,61106,61122,61140,61152,61157,61170,61198,
8 61217,61255,61274,61301,61313,61329,61357,61373,61390,61408,
9 61405,61423,61441,61443,61451,61465,61502,61510,61511,61515/
data (iy(i),i= 1601, 1700)/
+ 61526,61535,61553,61567,61573,61579,61599,61625,61632,61646,
1 61659,61671,61677,61687,61714,61723,61741,61774,61799,61825,
2 61834,61846,61854,61865,61875,61881,61894,61910,61928,61940,
3 61954,61969,61985,61990,62007,62027,62035,62047,62054,62064,
4 62073,62094,62121,62148,62160,62165,62170,62186,62199,62207,
5 62221,62247,62300,62319,62340,62357,62367,62382,62407,62422,
6 62430,62447,62471,62491,62518,62541,62564,62585,62613,62648,
7 62668,62697,62744,62768,62789,62809,62841,62869,62897,62920,
8 62950,62972,63003,63019,63046,63057,63067,63076,63089,63102,
9 63136,63171,63192,63220,63253,63276,63304,63328,63342,63350/
data (iy(i),i= 1701, 1800)/
+ 63359,63367,63381,63394,63408,63424,63438,63447,63459,63468,
1 63476,63490,63500,63566,63948,63946,63951,63957,63972,63989,
2 64032,64027,64025,64024,64028,64039,64058,64073,64119,64159,
3 64175,64203,64224,64236,64259,64288,64298,64309,64328,64349,
4 64364,64401,64426,64432,64440,64458,64483,64500,64521,64548,
5 64572,64589,64687,64788,64782,64785,64867,64910,65199,65372,
6 65447,65652,65819,65847,66090,66502,66333,66592,66615,66708,
7 67211,67556,67626,67732,67797,67909,68439,68457,68649,68664,
8 68679,68708,68758,68910,68960,68897,68909,69385,70177,70541,
9 70684,70590,71493,71937,72347,72949,72987,72990,74116,74121/
data (iy(i),i= 1801, 1900)/
+ 74415,74820,74875,74877,75859,75869,76960,77054,78320,79240,
1 80039,81231,82219,82618,82672,82726,83839,83994,84424,85049,
2 85824,85848,85903,85920,86263,87720,88459,88829,90009,90442,
3 90771,90786,90946,91446,91876,92225,92220,92792,93240,93497,
4 94039,94071,94134,94404,94430,94754,94814,94840,95043,95058,
5 95148,95241,95764,96003,96514,96520,96524,96529,96576,96589,
6 96864,97145,97152,97152,97152,97181,98022,97423,97100,97084,
7 96860,96811,96655,96516,96308,95960,95785,95673,95669,95356,
8 95151,95028,95013,95151,93611,93515,93413,93379,93321,92541,
9 92479,92049,91794,91798,91628,91625,91332,91153,90362,90344/
data (iy(i),i= 1901, 1920)/
+ 90340,90359,90371,90371,90940,91135,93166,96692,97808,98936,
1 99625,99703,99737,99737,99743,99737,99737,99737,99739,99740/
c
c map of british columbia
c this data was provided by jeff ovall
c
nvf=1920
nbf=nvf
ntf=3
rp(12)=0.1e0_rknd
rp(13)=1.5e0_rknd
c
sp(1)='British Columbia'
sp(2)='British Columbia'
sp(3)='British Columbia'
c
do i=1,nvf-1
if(ix(i) /= ix(i+1)) cycle
if(iy(i) /= iy(i+1)) cycle
write(6,*) 'error', i,ix(i),iy(i)
enddo
do i=1,nvf
vx(i)=real(ix(i),rknd)
vy(i)=real(iy(i),rknd)
enddo
c
c Vancouver Island
c
do i=1,285
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)=1
enddo
ibndry(2,285)=1
c
c Graham Island
c
do i=286,391
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)=2
enddo
ibndry(2,391)=286
c
c Mainland
c
do i=392,nbf
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)=3
enddo
ibndry(2,nbf)=392
c
eps=1.0e-3_rknd
call fixdg2(nvf,nbf,vx,vy,sf,ibndry,eps)
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 13.0 - - - september, 2018
c
c-----------------------------------------------------------------------
subroutine gd23(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), dimension(600), save :: ix,iy
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(i),i= 1, 100)/
+ 34609,35071,38054,38231,38260,38405,38494,39439,40660,40845,
1 40951,41434,42059,42748,42911,43039,43075,43728,43827,43991,
2 44239,44232,44012,44076,44140,44282,44438,44488,44516,44566,
3 44829,44758,45163,45191,45219,45305,45397,45468,45503,45511,
4 45511,45553,45788,45986,45993,46050,46057,46079,46086,46093,
5 46100,46128,46135,46185,46192,46199,46221,46228,46235,46228,
6 46235,46952,47172,47513,47641,47464,47527,47712,47883,47975,
7 48223,49495,49502,49509,49502,49097,49211,49225,49331,49353,
8 49360,49445,49480,49630,49878,49992,50070,50105,50091,50091,
9 50355,51225,51897,52031,52936,53025,53172,53468,54414,55071/
data (ix(i),i= 101, 200)/
+ 55398,55603,55905,56335,56576,56961,57109,57759,58461,58645,
1 59081,59207,59417,59528,59978,60411,60471,60489,60505,60537,
2 60563,61050,61178,61229,61390,61553,61680,61659,61624,61567,
3 61579,61797,61946,62186,62431,62491,62671,63276,63375,63469,
4 63565,63685,63774,64426,64814,65209,65472,65852,66065,66224,
5 66865,67254,67767,67954,68205,68393,68472,68591,68810,69099,
6 69199,70606,70739,70774,70682,70414,69906,69688,69609,69669,
7 69693,69732,69694,69401,69590,69968,70133,70218,70389,70422,
8 70390,70371,70275,70046,69965,70015,70226,70576,70771,71310,
9 71250,70841,70793,70824,70851,71140,71169,71134,70997,70709/
data (ix(i),i= 201, 300)/
+ 70199,70091,70030,70160,70085,70173,70218,70437,70882,70969,
1 70971,70854,70787,70756,70293,69918,69503,68876,68077,68054,
2 67892,67432,67334,67069,66941,66843,66738,66156,65787,65382,
3 65127,64945,64312,64006,63286,63120,63166,62878,62890,62827,
4 62601,62232,62198,62262,62240,62299,62488,62529,62538,62417,
5 62400,62663,62638,62381,62337,62219,62070,62016,62032,62204,
6 62148,61890,61753,61874,62039,62202,62393,62443,62411,62259,
7 62096,62074,62114,62325,62377,62314,62189,62102,62166,61865,
8 61939,61895,61704,61921,61751,61120,60367,60045,59464,59081,
9 58761,58589,58597,58415,58267,58299,58218,57803,57390,57173/
data (ix(i),i= 301, 400)/
+ 56681,56617,55431,55134,54931,54528,54136,53680,53326,53233,
1 53079,52943,52801,52474,52384,52078,51575,50938,50174,49910,
2 49715,49563,49478,49365,49246,49023,48971,48918,48830,48739,
3 48074,47904,47633,47605,47725,47640,47586,47424,47073,46907,
4 46610,46566,46578,46926,47065,47258,47289,47256,47208,47000,
5 46986,46971,47051,47318,47353,47403,47612,47697,47637,47432,
6 47438,47487,47403,47403,47486,47403,47136,47199,47403,47590,
7 48174,48298,48071,47403,46795,46480,46693,47187,47403,47894,
8 48545,49025,49374,49512,50075,50589,51693,52880,53345,53268,
9 53001,52966,51678,50588,49483,49252,48567,48053,47629,47403/
data (ix(i),i= 401, 500)/
+ 46849,46616,46501,46220,45996,46106,46011,45500,44933,44436,
1 44505,44793,45203,45838,46999,47403,48071,49115,50590,51130,
2 51485,52038,52516,53131,53507,54252,54889,55717,56565,56813,
3 57470,57873,59562,59559,59554,59551,59549,59458,59037,58070,
4 57775,57500,57369,57661,57912,58246,58802,59005,59133,59158,
5 59099,58986,58827,58164,57849,57746,57805,57939,58015,58013,
6 57873,57570,57495,57070,56629,56282,56085,55926,55918,56046,
7 55916,55427,55331,55045,54820,54819,54818,54650,54457,53830,
8 53694,53326,53118,53047,52876,52757,52690,52560,52520,52498,
9 51920,51744,51442,51391,51290,51363,51580,52175,52309,52228/
data (ix(i),i= 501, 572)/
+ 51659,51497,51584,51406,51097,51519,51898,49920,49753,49387,
1 49468,49176,48714,48687,48600,48657,48674,48600,48164,48164,
2 48231,48238,47903,47446,47403,47057,46782,46610,46542,46377,
3 46256,44475,44138,44016,43914,44155,44274,44159,43701,43718,
4 41103,41140,41156,40327,38126,37605,35821,35756,36284,35175,
5 35249,34509,32257,30798,28689,29839,30033,30977,31260,31254,
6 31190,31313,31413,32265,32469,33405,34831,35576,35803,35799,
7 35452,34351/
data (iy(i),i= 1, 100)/
+ 53333,54036,58922,58709,59029,59000,59299,61252,63943,64462,
1 64675,65783,67522,69390,69852,70306,70505,72863,73239,73978,
2 74745,75349,75349,75384,75377,75384,75384,75903,75974,76250,
3 77046,77174,78772,78885,79148,79439,79951,80306,80540,80703,
4 81094,81257,81414,82003,82067,82251,82308,82429,82500,82592,
5 82678,82990,83047,83359,83430,83509,83693,83750,83807,84275,
6 84332,84730,84595,84460,84453,84439,84346,84183,84112,84133,
7 84034,86413,86498,86598,86669,86882,87194,87322,87833,87983,
8 88082,88402,88558,89133,89978,90347,90752,90937,91058,91186,
9 91204,91214,91262,91272,91139,91172,91306,91386,91413,91490/
data (iy(i),i= 101, 200)/
+ 91529,91488,91429,90740,90508,90315,90279,90523,90570,90652,
1 91007,91152,91281,91291,91200,91337,91398,91565,91721,91780,
2 91829,92003,92259,92817,93438,93678,93963,94273,94407,94622,
3 94821,95239,95343,95609,96062,96117,95993,95268,95242,95318,
4 95520,95653,95688,95945,96140,96339,96552,97062,97224,97275,
5 97230,97286,97465,97598,97765,97961,98163,98814,99108,99285,
6 99345,100000,99914,99731,99406,99034,98548,98073,97610,97108,
7 96900,96573,96399,95977,95786,95516,95322,95198,94678,94377,
8 94249,94174,93966,93757,93620,93343,93279,93248,93171,92516,
9 91846,91541,91434,91359,91295,91019,90909,90846,90778,90749/
data (iy(i),i= 201, 300)/
+ 90863,90790,90656,90208,89403,89162,89038,87336,86815,86636,
1 86496,86223,86125,86080,85398,84701,84326,83376,82592,82677,
2 82755,82589,82560,82534,82442,82251,81903,81837,81522,81278,
3 81213,81112,81074,80908,80642,80482,80122,80025,79869,79807,
4 79754,79517,79254,79191,78795,78677,78589,78482,78430,78191,
5 77951,77490,77437,77440,77199,77201,77284,77259,77107,76871,
6 76790,76744,76581,76504,76508,76192,76126,76057,75927,75679,
7 75445,75315,75200,75113,75026,74901,74897,74830,74492,74412,
8 73835,73740,73693,73301,73315,73463,73761,73844,73866,73771,
9 73813,73943,74072,74271,74641,75152,75700,76295,76557,76602/
data (iy(i),i= 301, 400)/
+ 76609,76599,76414,76416,76462,76629,76874,77274,77438,77425,
1 77405,77348,77211,76898,76742,76491,76180,75906,75661,75576,
2 75452,75209,75074,74718,74063,73540,73221,72891,72819,72744,
3 72491,72309,71876,71187,70455,70265,70233,70134,69702,69385,
4 68813,68515,68218,67775,67528,66799,66471,66361,66203,65869,
5 65767,65665,65495,65253,65195,65114,64776,64482,64105,63724,
6 63552,63414,63356,63270,63084,62932,62445,61866,61631,61416,
7 61019,60609,60012,59910,59817,59419,58734,58554,58527,58467,
8 58666,59111,59287,59357,58913,58702,58322,58328,57912,56987,
9 56630,56583,56554,56174,55239,55043,54768,54551,54234,54025/
data (iy(i),i= 401, 500)/
+ 53513,53181,53015,52675,51715,50673,50227,49395,48767,47912,
1 47402,47005,46713,46567,46698,46784,46927,46852,46771,47011,
2 46941,47023,47058,47157,47350,48006,48620,48935,49467,49654,
3 50151,50273,50252,49809,49215,48834,48501,47986,46956,45465,
4 44950,44535,44019,43474,43183,42797,42154,41779,41246,40847,
5 40620,40191,39723,38828,38424,37903,37487,37170,36989,36461,
6 36350,36111,36068,35828,35578,34944,34568,34267,33804,33382,
7 32867,32715,32615,32317,31784,31294,30741,30403,30015,29025,
8 28456,28161,27824,27706,27419,27221,27109,26892,26616,26470,
9 24149,23691,23315,23252,22906,22791,22449,21464,21065,20760/
data (iy(i),i= 501, 572)/
+ 20215,19775,19265,18923,18327,16705,15249,12780,12203,10939,
1 8939, 8430, 7625, 7508, 7144, 6853, 6762, 6246, 4662, 3800,
2 3260, 2762, 2404, 2110, 2067, 1717, 1235, 807, 637, 631,
3 688, 0, 556, 984, 1729, 2427, 3126, 3606, 4303, 5535,
4 11933,12203,12320,14311,19721,21419,22432,23217,24509,25635,
5 27031,29694,33779,37462,43486,44228,44371,45072,45516,46205,
6 46942,47296,47583,48523,48692,49465,50653,50805,51151,51642,
7 52083,52958/
c
c map of israel
c this data was provided by jeff ovall
c
sp(1)='Israel'
sp(2)='Israel'
sp(3)='Israel'
c
nvf=572
nbf=nvf
ntf=1
rp(12)=0.1e0_rknd
rp(13)=1.5e0_rknd
c
do i=1,nvf
vx(i)=real(ix(i),rknd)
vy(i)=real(iy(i),rknd)
enddo
c
do i=1,nbf
ibndry(1,i)=i
ibndry(2,i)=i+1
ibndry(3,i)=0
ibndry(4,i)=2
ibndry(5,i)=0
ibndry(6,i)=1
enddo
ibndry(2,nbf)=1
c
eps=1.0e-3_rknd
call fixdg2(nvf,nbf,vx,vy,sf,ibndry,eps)
ip(1)=ntf
ip(2)=nvf
ip(3)=nbf
c
c make itnode
c
call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
return
end
.