p3.for

关于simple算法

*=======================================================================
      SUBROUTINE USER(K)
*-----------------------------------------------------------------------
$INCLUDE:'SIMPLE.INC'
      DIMENSION TH(22),THU(22),THDIF(22),THCV(22),THCVS(22)
      EQUIVALENCE(X,TH),(XU,THU),(XDIF,THDIF),(XCV,THCV),
     + (XCVS,THCVS),(XL,THL)
*-----------------------------------------------------------------------
      COMMON/CNTL/LSTOP
      COMMON/SORC/SMAX,SSUM
      COMMON/COEF/FLOW,DIFF,ACOF
*-----------------------------------------------------------------------
*---------------------------   PROBLEM  3   ----------------------------
*-----------------------------------------------------------------------
	DIMENSION T(22,22),W(22,22)
	EQUIVALENCE (F(1,1,4),W(1,1)),(F(1,1,5),T(1,1))
*-----------------------------------------------------------------------
	GOTO (1,2,3,4,5,6),K
*-----------------------------------------------------------------------
* 	ENTRY GRID
1	XL=0.5
	YL=0.5
	L1=7
	M1=7
	MODE=1
	CALL SUPPLY(1)
	RETURN
*-----------------------------------------------------------------------
*	ENTRY START
2	CONTINUE
	LSOLVE(4)=.TRUE.
	LPRINT(4)=.TRUE.
	LPRINT(5)=.TRUE.
	TITLE(4)=' W/WBAR'
	TITLE(5)=' THETA '
	LAST=22
	TW=0.
	DO 100 J=1,M1
	DO 100 I=1,L1
	W(I,J)=0.
	T(I,J)=1.
	T(L1,J)=TW
	T(I,M1)=TW
100	CONTINUE
	AMU=1.
	DEN=1.
	COND=1.
	CP=1.
	RHOCP=DEN*CP
	DPDZ=-100.
	DTBDZ=5.
	RETURN
*-----------------------------------------------------------------------
*	ENTRY DENSE
3	CONTINUE
	RETURN
*-----------------------------------------------------------------------
*	ENTRY BOUND
4	CONTINUE
	ASUM=0.
	WSUM=0.
	TSUM=0.
	DO 300 J=2,M2
	DO 300 I=2,L2
	AR=XCV(I)*YCV(J)
	WR=W(I,J)*AR
	ASUM=ASUM+AR
	WSUM=WSUM+WR
300	TSUM=TSUM+WR*T(I,J)
	WBAR=WSUM/ASUM
	TB=TSUM/(WSUM+1.E-30)
	DH=4.*XL*YL/(XL+YL)
	RE=DEN*WBAR*DH/AMU
	FRE=(-DPDZ*DH/(0.5*DEN*WBAR**2+1.E-30))*RE
	QW=DTBDZ*RHOCP*WSUM/(XL+YL)
	ANU=QW*DH/(COND*(TW-TB)+1.E-30)
	IF(ITER.NE.10)GOTO 301
	LSOLVE(4)=.FALSE.
	LSOLVE(5)=.TRUE.
301	CONTINUE
	RETURN
*-----------------------------------------------------------------------
*	ENTRY OUTPUT
5	CONTINUE
	IF(ITER.NE.0)GOTO 400
	WRITE(*,401)
401	FORMAT('   ITER',12X,'F*RE',17X,'NU')
400	WRITE(*,402)ITER,FRE,ANU
402	FORMAT(I6,1P2E20.4)
	IF(ITER.NE.LAST)RETURN
	DO 410 J=1,M1
	DO 410 I=1,L1
	W(I,J)=W(I,J)/WBAR
410	T(I,J)=(T(I,J)-TW)/(TB-TW)
	CALL SUPPLY(2)
	RETURN
*-----------------------------------------------------------------------
*	ENTRY GAMSOR
6	CONTINUE
	DO 500 I=1,L1
	DO 500 J=1,M1
	GAM(I,J)=AMU
	IF(NF.EQ.5)GAM(I,J)=COND
	GAM(I,1)=0.
	GAM(1,J)=0.
500	CONTINUE
	IF(NF.EQ.5)GOTO 511
	DO 510 J=2,M2
	DO 510 I=2,L2
	CON(I,J)=-DPDZ
510	CONTINUE
	RETURN
511	DO 520 J=2,M2
	DO 520 I=2,L2
	THETA=(T(I,J)-TW)/(TB-TW+1.E-30)
	DTDZ=THETA*DTBDZ
520	CON(I,J)=-RHOCP*W(I,J)*DTDZ
	RETURN
	END