浮遊砂の河床変動計算できた

 今夜3日目にして、浮遊砂の計算モデルが完成しました。以前挫折してホームページに中途半端なソースを載せたのとは、違う形になりました。初めに初期水面形を計算し、2回目以降ルーチンの計算を行う。そして最後に出力を行ってプログラムは終了。
 不等流計算で、差を一回目計算してから二回目以降のルーチンを計算するようにしたのと同じ感じだ。これで何とか現場のための水理学を終えることができる。下にプログラム本体を載せる。サブルーチンは、資料に載っているものとほぼ同じだ。ただ、wfの計算だけはうまくいかなかったので、計算値を定数として使った。
 課題は、流砂量式の選び方と差分式の導き方だ。


program main
! use globals
use huyuu
implicit none
integer:: etim=86400,j
real(8):: ramda=0.4,cc(60)

!------------------------------------------------------------------------
open(20,file='input z.csv')
do j=1,nj
read(20,*)z(j)
enddo
close(20)
!-----------------------------------
ot(1)=0
ot(2)=1
ot(3)=4
ot(4)=8
ot(5)=12
ot(6)=24
!------------------------------------

do j=1,6
ptim(j)=ot(j)*3600
enddo
do j=1,6
! write(*,*)ptim(j)
enddo
!------------------------------------
do j=1,nj
x(j)=200.*(j-1)
z(j)=200./700*(j-1)
enddo
do j=1,nj
zz(j)=z(j)
enddo

!*************************************************************************


tim2=0
itm=1
call caliwa(s,g,d,usc,tsc)
!-------------------------------------------

call calh(h,z,b,q,sn,h0,g,dx,eps,ie,nj)
do j=1,nj
ts(j)=h(j)*ie(j)/(s*d)
us(j)=sqrt(g*h(j)*ie(j))
enddo
!
dz(nj)=0.

call calqsu(s,g,d,nj,us,ts,tsc,wf,qsu)
do j=1,nj-1
cc(j)=0.
enddo
call calc(nj,x,b,qsu,wf,c,q,h,idt1,idt2,tim,cc)
call pout(z,zz,h,ts,us,qsu,c,nj,tim2,ptim,itm)
call caldz(nj,ramda,idt2,wf,z,qsu,c,dz)

write(*,*)tim2
1000 continue
!
open(30,file='out.d')
do j=1,nj
! write(30,*)z(j)
enddo
close(30)
open(40,file='out.d')
do j=1,nj
! read(40,*)z(j)
enddo
close(40)

write(*,*)'calh start'

call calh(h,z,b,q,sn,h0,g,dx,eps,ie,nj)
write(*,*)'calh end'
do j=1,nj
ts(j)=h(j)*ie(j)/(s*d)
us(j)=sqrt(g*h(j)*ie(j))
enddo
!
call calqsu(s,g,d,nj,us,ts,tsc,wf,qsu)
write(*,*)'calqsu'
call calc(nj,x,b,qsu,wf,c,q,h,idt1,idt2,tim,cc)
write(*,*)'calc'
call pout(z,zz,h,ts,us,qsu,c,nj,tim2,ptim,itm)
call caldz(nj,ramda,idt2,wf,z,qsu,c,dz)
write(*,*)'caldz'
!
tim2=tim2+idt2
write(*,*)tim2
if(tim2.lt.etim) goto 1000
call pout(z,zz,h,ts,us,qsu,c,nj,tim2,ptim,itm)

end program main




ブログ気持玉

クリックして気持ちを伝えよう!

ログインしてクリックすれば、自分のブログへのリンクが付きます。

→ログインへ

なるほど(納得、参考になった、ヘー)
驚いた
面白い
ナイス
ガッツ(がんばれ!)
かわいい

気持玉数 : 1

なるほど(納得、参考になった、ヘー)

この記事へのコメント

この記事へのトラックバック