关于un上一步长数值的问题
本帖最后由 独孤 于 2009-2-23 14:29 编辑我在求解温度场，其中材料参数是变物性的，需要根据上一个时间步长的温度来决定下一个时间步长的材料参数，下面是我的PDE文件。
我发现程序求解需要的时间是不变材料参数时候程序运行时间的N倍。这个N是PRE里面定义的材料数目。我定义了51个材料（包括15个边界），这就意味着程序运算量一下子增长了50倍。不明白这是为什么？
disp u
coor x y z
coef un
mate rou ec ek eq wb tb hz eqb 4.0d-6;1.04d-5;0.0d0;4.0d-6;1;37;0;0;
shap %1 %2
gaus %3
mass %1 rou*ec
stif
$c6 if (imate.gt.1) then
$cv eqb=un
$c6 endif
$c6 if (fxy.gt.50) then
$cv eqb=37
$c6 endif
$c6 if (fxy.lt.0) then
$cv eqb=37
$c6 endif
dist =
+*ek
+*ek
+*ek
load =+【u】*(eq+4077.15*wb*(tb-eqb)+hz)
end
[ 本帖最后由 独孤 于 2008-12-29 14:23 编辑 ] 原因就在于你插入了这些行的程序，因为这些程序是对每一个单元要进行循环的，你看一下e元件程序以及其调用的子程序就可以了。
$c6 if (imate.gt.1) then
$cv eqb=un
$c6 endif
$c6 if (fxy.gt.50) then
$cv eqb=37
$c6 endif
$c6 if (fxy.lt.0) then
$cv eqb=37
$c6 endif
再一个你的pde里面的fxy找不到出处，好像不对，具体应该imate吧。 因为“插入了这些行的程序”？？如果我不插入这些$c6的程序行，就不用对每一个单元循环？我把PDE修改成这样，
disp u
coor x y z
coef un
mate rou ec ek eq wb tb hz eqb 4.0d-6;1.04d-5;0.0d0;4.0d-6;1;37;0;0;
shap %1 %2
gaus %3
mass %1 rou*ec
stif
dist =
+*ek
+*ek
+*ek
load =+*(eq+4077.15*wb*(tb-un)+hz)
end
我对程序运行，发现还是要循环51次（我定义了51个材料）
它的程序我看得不是太明白，有一点疑问：
我仅仅是在方程+*(eq+4077.15*wb*(tb-un)+hz)中调用了UN，也就是上一个时间步长的单元的温度而已，为什么还是要循环51次？ 下面是e元件程序中的一部分，你看一下就明白了，呵呵。
numtyp = 3 !你的肯定是51
IF (IT.EQ.0) THEN
open(11,file='unod',form='unformatted',status='old')
ELSE
open(11,file='unod',form='unformatted',status='old')
ENDIF
read(11) ((eu1(j,i),i=1,knode),j=1,kdgof)
C.......OPEN DIAG file
OPEN (2,FILE=' ',FORM='UNFORMATTED',STATUS='OLD')
READ(2) (numcol(I),I=1,NEQ1)
READ(2) (na(I),I=1,maxa)
CLOSE(2)
C.......OPEN ELEM0 file
OPEN (3,FILE=' ',FORM='UNFORMATTED',STATUS='OLD')
itime=0
1 continue
itime=itime+1
if (itime.gt.1) then
write(*,*) 'Nonlinear Iteration Times ========',itime
rewind(3)
rewind(32)
endif
READ (32) ((U(J,I),J=1,KDGOF),I=1,KNODE)
cc WRITE (*,*) 'BF ='
cc WRITE(*,7) ((U(J,I),J=1,KDGOF),I=1,KNODE)
DO 111 I=1,MAXA
A(I) = 0.0
111 CONTINUE
DO 2300 I=1,NEQ
2300 CONTINUE
NUMEL=0
DO 2000 ITYP=1,NUMTYP ！这里要对材料进行循环求解 51次，所有有51个循环了
C .... READ ELEMENT NODAL CONNECTION POINTS AND MATERIAL No.
READ (3) NUM,NNODE,
* ((NODE((I-1)*NNODE+J),J=1,NNODE),I=1,NUM)
cc WRITE(*,*) 'NUM =',NUM,' NNODE =',NNODE
cc WRITE(*,*) 'NODE ='
cc WRITE(*,6) ((NODE((I-1)*NNODE+J),J=1,NNODE),I=1,NUM)
NNE = NNODE
nne = nne-1
K=0
DO 115 J=1,NNE
JNOD = NODE(J)
DO 115 L=1,KDGOF
IF (NODVAR(L,JNOD).NE.0) K=K+1
115 CONTINUE
WRITE(*,*) 'K =',K
........
为节省篇幅，以下略
版主太强悍了！ 啥都不说了，Happy New Year! 不知道有没有办法不让程序进行51次循环？能否通过修改版主提出的那段代码来实现？
祝大家牛年吉祥！ 好奇怪，我的numtyp = 2，不是51
这是为什么？
原帖由 zhoulon 于 2008-12-31 12:44 发表 http://forum.simwe.com/images/common/back.gif
下面是e元件程序中的一部分，你看一下就明白了，呵呵。
numtyp = 3 !你的肯定是51
IF (IT.EQ.0) THEN
open(11,file='unod',form='unformatted',status='old')
ELSE
... numtyp是单元类型数，＝2是对的
原帖由 独孤 于 2009-2-13 13:16 发表 http://forum.simwe.com/images/common/back.gif
好奇怪，我的numtyp = 2，不是51
这是为什么？
请问能否实现不让程序因为材料数目的缘故而迭代循环N次？
原帖由 独孤 于 2009-2-13 13:16 发表 http://forum.simwe.com/images/common/back.gif
好奇怪，我的numtyp = 2，不是51
这是为什么？
我上面说错了，numtyp是单元类型数目，你的程序应该有一种体单元，一种面单元，所以这个值是2了。
具体参数在e元件程序里面的prmt数组里面，你自己看一下吧。 除了fxy 应该是imate之外，可能别的东东有问题，非线性问题收敛的迭代步数跟算法有关
原帖由 独孤 于 2009-2-13 13:48 发表 http://forum.simwe.com/images/common/back.gif
请问能否实现不让程序因为材料数目的缘故而迭代循环N次？
本帖最后由 独孤 于 2009-2-14 08:33 编辑
seeksilence 说的是对的。从bat文件的循环上面能够发现：
Par.dat文件内容：
starta disp0 id0 nv disp0 bfd time0 time coor0 unod para.io elem0 sys diag
if exist stop del stop
:1
bft disp0 time coor0 nv bfdt stop
IF EXIST END DEL END
:2
epara disp0 sys time nv coor0 bfd diag elem0 f displace end
IF NOT EXIST END GOTO 2 ----- 这个判断语句是进行N次迭代（在此例中是51次）的原因
call post.bat
if not exist stop goto 1
Par.gcn文件内容如下：
defi
a npar &
STARTsin a
if exist stop del stop
:1
bft
IF EXIST END DEL END
:2
SOLVsin a
IF NOT EXIST END GOTO 2
post a
if not exist stop goto 1
请问大侠，是不是我的算法选择有问题？
PS：
par.pde 文件的内容经过修改为：
disp u
coor x y z
coef un
mate rc ek eq wb tb hz 1085;0.42;600;0.001;37;0;
shap %1 %2
gaus %3
mass %1 rc
stif
dist =
+*ek
+*ek
+*ek
load =+*(eq+3901.85*wb*(tb-un)+hz)
end 本案例求解对象是温度场，没有任何耦合，仅仅存在非线性的特点。请问能否推荐下运算速度快的求解器？
重要的是，我能否通过修改.gcn文件中的求解器来一个一个地做试验选定最好的求解器？不知道这个思路是否可行！
祝情人节快乐！ 给的材料参数是否合理？ （rc ek 1085;0.42;）
先别管快不快，分析下结果合不合理，结果不对再快也没意义 有个getsolv命令可以修改求解器的。
页:
[1]