c*********************** problem name: battery ************************
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
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)
integer(kind=iknd), save, dimension(5) :: ic
real(kind=rknd), dimension(*) :: values
common /val0/k0,ku,kx,ky,kl,kuu,kux,kuy,kul,
+ kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll
cy
data ic/2500,700,500,20,5/
c
coeff=real(ic(itag),rknd)*1.0e-2_rknd
values(k0)=coeff*ux
values(kx)=coeff
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
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)
integer(kind=iknd), save, dimension(5) :: ic
real(kind=rknd), dimension(*) :: values
common /val0/k0,ku,kx,ky,kl,kuu,kux,kuy,kul,
+ kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll
cy
data ic/250000,8000,1,2000,500/
c
coeff=real(ic(itag),rknd)*1.0e-4_rknd
values(k0)=coeff*uy
values(ky)=coeff
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
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)
integer(kind=iknd), save, dimension(5) :: ic
real(kind=rknd), dimension(*) :: values
common /val0/k0,ku,kx,ky,kl,kuu,kux,kuy,kul,
+ kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll
cy
data ic/0,-1,-1,0,0/
c
c
values(k0)=real(ic(itag),rknd)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
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
real(kind=rknd), save, dimension(4) :: g,al
common /val1/k0,ku,kl,kuu,kul,klu,kll
cy
data g/0.0e0_rknd,1.0e0_rknd,2.0e0_rknd,3.0e0_rknd/
data al/0.0e0_rknd,3.0e0_rknd,2.0e0_rknd,1.0e0_rknd/
c
values(k0)=g(itag)-al(itag)*u
values(ku)=-al(itag)
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
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,kll,klb,kub,kic,kim,kil
cy
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
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,kuu,kux,kuy,kul,
+ kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll
cy
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
c
c-----------------------------------------------------------------------
subroutine p2xy(x,y,dx,dy,u,ux,uy,rl,itag,jtag,values)
c
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,kuu,kux,kuy,kul,
+ kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll
cy
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
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
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
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(*) :: values
common /val4/jx,jy,jxs,jys,jxl,jyl,jxss,jyss,jxll,jyll,
+ jxsl,jysl,jxls,jyls
cy
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
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
return
end
c-----------------------------------------------------------------------
c
c piecewise lagrange triangle multi grid package
c
c edition 11.0 - - - june, 2012
c
c-----------------------------------------------------------------------
subroutine gdata(vx,vy,sf,itnode,ibndry,ip,rp,sp,iu,ru,su,sxy)
cx
use mthdef
implicit real(kind=rknd) (a-h,o-z)
implicit integer(kind=iknd) (i-n)
integer(kind=iknd), dimension(5,*) :: itnode
integer(kind=iknd), dimension(7,*) :: ibndry
integer(kind=iknd), dimension(100) :: ip,iu
integer(kind=iknd), save, dimension(32) :: label,it1,it2
integer(kind=iknd), save, dimension(5) :: ix
integer(kind=iknd), save, dimension(9) :: iy
integer(kind=iknd), save :: ntf,nvf,nbf,nx,ny,ifirst
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
c
c thermal battery problem
c data courtesy of l. demcowicz
c
data ifirst/1/
data nx,ny/5,9/
data ix/0,61,65,80,84/
data iy/0,8,16,36,120,188,212,232,240/
data ntf,nvf,nbf/32,45,76/
data hmax,grade/0.1e0_rknd,1.5e0_rknd/
data label/2,3,3,2,5,5,1,1,1,1,4,5,4,4,4,4,
+ 1,5,5,5,5,5,5,1,1,1,1,1,1,1,1,1/
data it1/16,16,21,26,31,11,36,6,2,42,7,37,12,17,22,27,
+ 8,8,13,18,23,28,33,38,4,44,9,14,19,24,29,34/
data it2/11,44,45,46,47,6,48,1,50,37,51,32,52,53,54,55,
+ 3,60,61,62,63,64,65,66,68,39,69,70,71,72,73,74/
c
c 1 1 1 1 8 10 24 28
c 5 5 5 1 6 12 23 32
c 2 4 5 1 1 16 22 31
c 3 4 5 1 2 15 21 30
c 3 4 5 1 3 14 20 29
c 2 4 5 1 4 13 19 28
c 5 4 5 1 5 11 18 27
c 1 1 1 1 7 9 17 25
c
if(ip(41)==1) then
sp(2)='battery'
sp(1)='battery'
sp(3)='battery'
sp(6)='battery_mpixxx.rw'
sp(7)='battery.jnl'
sp(9)='battery_mpixxx.out'
endif
c
nvf=0
do j=1,ny
do i=1,nx
nvf=nvf+1
vx(nvf)=real(ix(i),rknd)*1.0e-1_rknd
vy(nvf)=real(iy(j),rknd)*1.0e-1_rknd
enddo
enddo
c
nbf=0
do j=1,ny-1
do i=1,nx
nn=(j-1)*nx+i
nbf=nbf+1
ibndry(1,nbf)=nn
ibndry(2,nbf)=nn+nx
ibndry(3,nbf)=0
ibndry(4,nbf)=0
if(i==1.or.i==nx) ibndry(4,nbf)=1
ibndry(5,nbf)=0
ibndry(6,nbf)=0
ibndry(7,nbf)=0
if(i==1) ibndry(7,nbf)=1
if(i==nx) ibndry(7,nbf)=3
enddo
enddo
c
do i=1,nx-1
do j=1,ny
nn=(j-1)*nx+i
nbf=nbf+1
ibndry(1,nbf)=nn
ibndry(2,nbf)=nn+1
ibndry(3,nbf)=0
ibndry(4,nbf)=0
if(j==1.or.j==ny) ibndry(4,nbf)=1
ibndry(5,nbf)=0
ibndry(7,nbf)=0
if(j==1) ibndry(7,nbf)=4
if(j==ny) ibndry(7,nbf)=2
enddo
enddo
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
ip(5)=max(ip(5),ifirst)
ip(6)=1
ip(19)=1
ip(20)=5
rp(15)=hmax
rp(16)=grade
c
c make itnode (saved output as it1/it2 because
c we get different ordering for single and double
c and label depends on the ordering)
c
cc call sklutl(0_iknd,vx,vy,sf,itnode,ibndry,ip,rp,iflag,sxy)
c
c label regions
c
do i=1,ntf
itnode(1,i)=it1(i)
itnode(2,i)=it2(i)
itnode(5,i)=label(i)
enddo
return
end
.