LS-DYNA自定义材料本构子程序的简便编译方法

gdyu_yu

df_af_aq曾经介绍过ls-dyna自定义本构子程序的基本操作流程，
（http://forum.simwe.com/thread-994915-1-1.html）
按这个流程做下来确实可以生成自定义求解器。但客观地说，这个流程略显繁琐。下面介绍一种简化编译方法，以供大家参考。
基本编译思想：写好fortran子程序之后，通过IVF编译器以及nmake命令编译for子程序，生成obj目标文件和求解器。
简化方法：
1、利用文本编辑器写fortran子程序，起个方便查看的名字, 存为.f格式文件，比如umat41-for-shell.f。
2、打开dyn21.f文件，注释掉subroutine umat41 (cm,eps,sig,epsp,hsv,dt1,capa,etype,tt,temper,failel,crv,cma)下这一整个子程序的语句（包括这个语句），或者直接删掉这部分子程序的内容。保存dyn21.f。
3、用文本编辑器打开makefile，做如下修改：
1）  OBJS= dyn21.obj dyn21b.obj一行中，增加umat41-for-shell.obj。变为OBJS= dyn21.obj dyn21b.obj umat41-for-shell.obj
2）参照dyn21.obj: dyn21.F
$(FC) $(FFLAGS) dyn21.F
格式，书写两行语句
umat41-for-shell.obj：umat41-for-shell.f
$(FC) $(FFLAGS) umat41-for-shell.f
比如：
dyn21.obj: dyn21.F
$(FC) $(FFLAGS) dyn21.F
dyn21b.obj: dyn21b.F
$(FC) $(FFLAGS) dyn21b.F
umat41-for-shell.obj: umat41-for-shell.f
$(FC) $(FFLAGS) umat41-for-shell.f
然后保存makefile。
3）用IVF编译器和nmake命令编译子程序。这步以及后续的步骤与df_af_aq写的基本流程一样。

简化编译方法的优点：不需要每次都进入dyn21.f中修改子程序，也可以按fortran语法命名所需编译子程序的名称；编译简单快捷。

至于版上有些朋友反映，说赵海欧先生一书那个子程序无法输出应力，这大多情况是因子程序的变量名与971所要求的不一致引起的。
赵先生一书所用的版本是970，而根据关键字用户手册附录A的说明，不同的求解器版本，其子程序的变量名很可能不同。
如果你用971求解器所需的lib文件包，那么在写子程序的时候，尤其需要注意其变量名。
比如，在971中，前一时间步的历史变量数组名为hsv，而不是hisv；971中增加了内置历史变量epsp（前一时间步的等效塑性应变）；
当前时间步长名为dt1，当前时间为tt。所有这些变量以及这些变量在子程序中的定义、引用都需要与971的要求保持一致。
建议想二次开发LS-DYNA材料模型子程序的朋友，多看看LS-DYNA keyword user manual中APPENDIX A:  User Defined Materials部分的介绍。

笨笨驴吃萝卜

      昨天用LS-dyna 971.r511版本运行了赵海鸥先生编著的那本书里的程序，终于有结果了，之前按照版主的方法编译连接成功，但一直没有结果，后来在他们的官网上看见
http://ftp.lstc.com/anonymous/ou ... fined_materials.faq

Question:
Is there a difference in UMAT41-UMAT50?
Is the routine name of UMAT only connected with the parameter MT in
*MAT_USER_DEFINED_MATERIAL_MODEL?

子程序41-50是否存在区别？子程序名是否只与关键字中，*MAT_USER_DEFINED_MATERIAL_MODEL的参数MT有关？
Answer:
Avoid using umat42 and the vectorized routines umat48v and umat49v; these are special, undocumented subroutines that do not reside in dyn21.f and thus the user cannot modify them.  Aside from those exceptions, it doesn't matter which MT (from 41 to 50) you choose.  
The umat subroutine used is controlled directly by MT and by IVECT.
If IVECT is set to 1, the "v" (vectorized) version of the umat subroutine will be called, e.g., subroutine umat41v is called if MT=41 and IVECT=1.
在用户自定义程序中避免使用umat42子程序及向量子程序umat48v和umat49v，因为这些子程序是特殊的。

    我当时在dyn21.f文件里搜索发现umat42程序已经被引用掉了，怕麻烦，所以就直接用了42号材料，结果才发现因为42号等三个子程序是特殊的，在自定义程序里面是不能改动，所以一直没有结果，下面是我按照ls-dyna 971关键字手册和材料自定义dyn21.f程序包的函数变量改的子程序，已经运行通过了。注意要把dyn21.f中响应编号的材料自定义程序给引用掉，要不然编译会提示有错误的。另外就是，一定在自编程序中要include两个库函数，不然编译也不会成功

zdsimwe

第一编译容易出问题，大家可以用个简单的方法来测试一下：
1）把下载的开发文件放到一个目录下，包括所有的Fortran源文件和库文件。
2）选Intel Parallel Studio XE ... -> Command Prompt -> Parallel Studio XE with Intel Compiler XE .. -> Intel 64 Visual Studio 2010 mode

在开发文件的目录下，运行nmake
测试的时候直接用41号材料，线弹性，或者45号材料，超弹性。

在什么都不改的情况下，这两个材料都能算，应力位移都正常。然后才开始自己的。。。。。

zdsimwe

是少MKL库，但是你的编译器版本不对，出的问题谁也不知道。R7以前的都是ifort 10.1，之后的都是ifort 13.1 不能混淆。

lojade

老大V5  87
最近也在准备朝自己开发程序方面发展，还没有下定决心。先看看老大的作品。

76243386

先赞一个，斑竹真牛，大哥的神贴！

heraculia

谢谢版主的回复，已阅读您的大作，受益匪浅。
刚从abaqus的材料子程序转到dyna这边来做，以后还要多多向版主请教~~~

802

可以麻烦斑竹大哥把fortran10.1和MVS2008传给我一份吗？是在下不到，头疼的很。。。邮箱：654900815@qq.com~非常感谢！

dark-f

这个以后具体实现时会用到，谢谢先！

笨笨驴吃萝卜

经过我的运行发现利用 visual studio 2005+ IVF11.1.038也可以编译通过。

gdyu_yu

笨笨驴吃萝卜 发表于 2013-10-26 18:56
经过我的运行发现利用 visual studio 2005+ IVF11.1.038也可以编译通过。

事实上，没有MVS，只有IVF好像也是可以编译的。用好makefile即可。

hongjian

斑竹您好！我有一个ls-dyna子程序的问题想请教您，希望您不吝指教。
子程序nmake产生的ls971和ls971_d_R5.1.1_winx64_p或者ls971_s_R5.1.1_winx64_p有区别吗？因为我算了一个例子，温度相关的，ls971_d_R5.1.1_winx64_p可以算通，但用ls971算不通（dyna21.F未做任何改动，k文件不含子程序）。
请问您有注意过这个情况吗？

gdyu_yu

hongjian 发表于 2013-12-17 20:59
斑竹您好！我有一个ls-dyna子程序的问题想请教您，希望您不吝指教。
子程序nmake产生的ls971和ls971_d_R5.1 ...

我没有比较过。可能与温度有关的问题，需要用双精度求解器吧。

惆怅依旧

笨笨驴吃萝卜 发表于 2013-10-28 10:23
昨天用LS-dyna 971.r511版本运行了赵海鸥先生编著的那本书里的程序，终于有结果了，之前按照版主的方 ...

我用你这个编译，提示有错误，能帮我看看错在哪吗？
C:\Documents and Settings\Administrator.NEU-B02BEECCC67\桌面

惆怅依旧

笨笨驴吃萝卜 发表于 2013-10-28 10:23
昨天用LS-dyna 971.r511版本运行了赵海鸥先生编著的那本书里的程序，终于有结果了，之前按照版主的方 ...

我用你这个编译，提示有错误，能帮我看看错在哪吗？
提示：

E:\>c Updating The Yield Stress
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      ak=cm(6)+qh*hsv(1)
'ak' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Begin strain rate effect
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      if(cm(8).ne.0) then
此时不应有 then。

E:\>      d1d=eps(1)+davg
'd1d' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      d2d=eps(2)+davg
'd2d' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      d3d=eps(3)+davg
'd3d' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      d4d=eps(4)
'd4d' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      d5d=eps(5)
'd5d' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      d6d=eps(6)
'd6d' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      ds=d4d*d4d+d5d*d5d+d6d*d6d
'ds' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Computing Effective Strain
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      effs=sqrt(2.*(d1d*d1d+d2d*d2d+d3d*d3d+2.*ds)/3.)
'effs' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Computing Strain Rate
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      if (tt.ne.0) then
此时不应有 then。

E:\>      effs=effs/dt1
'effs' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      endif
'endif' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Computing the Strain rate Sensitivity using Cowper-Symond
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      powp=1/cm(9)
'powp' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      bas=effs/cm(8)
'bas' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      mult1=bas**powp
'mult1' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      mult=1.+bas**powp
'mult' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      endif
'endif' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      ak=mult*ak
'ak' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c End strain rate effect
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Computing Deviatoric Stress
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      q1=hsv(2)+sig(1)"Documents and Settings"! 3
'q1' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      q2=hsv(2)+sig(2)
'q2' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      q3=hsv(2)+sig(3)
'q3' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      q4=sig(4)
'q4' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      q5=sig(5)
'q5' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      q6=sig(6)
'q6' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      aj2=q4*q4+q5*q5+q6*q6+(q1*q1+q2*q2+q3*q3)/2
'aj2' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Computing Yield Function
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      ak2=3*aj2-ak*ak! 4
'ak2' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      scle=0
'scle' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Checking Yield
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      if(ak2.gt.0) then! 5
此时不应有 then!。

E:\>      scle=1
'scle' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      endif
'endif' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      fac1=1/(3.0*g+qh)
'fac1' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      fac2=3.0*g
'fac2' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      aj1=sqrt(3*abs(aj2))+1-scle
'aj1' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Computing plastic strain Increment
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      depi=scle*fac1*(aj1-ak)config.sys! 6
'depi' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Computing Total Plastic Strain
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      hsv(1)=hsv(1)+depi
'hsv' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      deps=scle*fac2*depi/aj1
'deps' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>c Stress Update
'c' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      sig(1)=sig(1)-deps*q1! 7
'sig' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      sig(2)=sig(2)-deps*q2
'sig' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      sig(3)=sig(3)-deps*q3
'sig' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      sig(4)=sig(4)-deps*q4
'sig' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      sig(5)=sig(5)-deps*q5
'sig' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      sig(6)=sig(6)-deps*q6
'sig' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      endif
'endif' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      return
'return' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>      end
'end' 不是内部或外部命令，也不是可运行的程序
或批处理文件。

E:\>cd ls971win32

E:\ls971win32>nmake

Microsoft (R) Program Maintenance Utility   Version 6.00.8168.0
Copyright (C) Microsoft Corp 1988-1998. All rights reserved.

        ifort -c -W0 -WB -unroll -fp:precise -4Yportlib -assume:byterecl,buffere
d_io  -Qfpp2 -DPCWIN -DINTEL -QxW -DIA32 -DIA32ONLY -DLSTCODE -DLONGFORMAT -DOPE
NMP -Qopenmp -DDOUBLEXR dyn21.F
Intel(R) Visual Fortran Compiler XE for applications running on IA-32, Version 1
2.0.0.104 Build 20101006
Copyright (C) 1985-2010 Intel Corporation.  All rights reserved.
ifort: command line remark #10010: option '/unroll' is deprecated and will be re
moved in a future release. See '/help deprecated'
ifort: command line warning #10212: /fp:precise evaluates in source precision wi
th Fortran.
ifort: command line remark #10010: option '/Qfpp2' is deprecated and will be rem
oved in a future release. See '/help deprecated'
ifort: command line remark #10279: option '/QxW' is deprecated and will be remov
ed in a future release. See '/help deprecated'

dyn21.F(2044): error #6353: A RETURN statement is invalid in the main program.
      return
------^
dyn21.F(1941): error #6221: The assumed-size array must be a dummy argument.   [
CM]
      dimension cm(*),eps(*),sig(*),hsv(*),crv(lq1,2,*),cma(*)
----------------^
dyn21.F(1941): error #6221: The assumed-size array must be a dummy argument.   [
EPS]
      dimension cm(*),eps(*),sig(*),hsv(*),crv(lq1,2,*),cma(*)
----------------------^
dyn21.F(1941): error #6221: The assumed-size array must be a dummy argument.   [
SIG]
      dimension cm(*),eps(*),sig(*),hsv(*),crv(lq1,2,*),cma(*)
-----------------------------^
dyn21.F(1941): error #6221: The assumed-size array must be a dummy argument.   [
HSV]
      dimension cm(*),eps(*),sig(*),hsv(*),crv(lq1,2,*),cma(*)
------------------------------------^
dyn21.F(1941): error #6917: This array with adjustable declarators must be decla
red a dummy argument.   [CRV]
      dimension cm(*),eps(*),sig(*),hsv(*),crv(lq1,2,*),cma(*)
-------------------------------------------^
dyn21.F(1941): error #6221: The assumed-size array must be a dummy argument.   [
CMA]
      dimension cm(*),eps(*),sig(*),hsv(*),crv(lq1,2,*),cma(*)
--------------------------------------------------------^
compilation aborted for dyn21.F (code 1)
NMAKE : fatal error U1077: 'ifort' : return code '0x1'
Stop.

E:\ls971win32>

hongjian

gdyu_yu 发表于 2013-12-18 13:12
我没有比较过。可能与温度有关的问题，需要用双精度求解器吧。

这个小小的bug解决了。确实应该是和精度有关系，我之前下载的是r511双精度的，什么也不修改nmake产生的ls971也不能对温度的k文件进行计算。刚刚新下了一个单精度的r511，倒是算过去了。谢谢版主的提示！
具体原因真比较迷糊，不知道版主你有单精度和双精度区别这方面的知识或者资料吗？

cainiao07

难怪我试的赵海鸥的例子，算的进入死循环了。看来还真是971版变量名有变化。

笨笨驴吃萝卜

惆怅依旧 发表于 2013-12-18 19:40
我用你这个编译，提示有错误，能帮我看看错在哪吗？
提示：

我现在从你给的错误信息提示中看到给我的感觉是你在子程序中没有定义上述一些变量，所有子程序用到的变量和数组都需要在子程序里面定义。后面的就是你是不是没有把dyna21.F中给出相应的usemat41（具体不记得了，你找一下）子程序是不是没有注释掉，造成重复定义了。

笨笨驴吃萝卜

惆怅依旧 发表于 2013-12-18 19:40
我用你这个编译，提示有错误，能帮我看看错在哪吗？
提示：

赵海鸥先生书里面给的子程序名称一定要和dyn21.F给出的子程序名称不能有重复定义，不然就会有上面的提示。

惆怅依旧

笨笨驴吃萝卜 发表于 2013-12-25 20:37
我现在从你给的错误信息提示中看到给我的感觉是你在子程序中没有定义上述一些变量，所有子程序用到的变量 ...

嗯，我没有删完，因此就出现错误了

guoshuai0501

楼主厉害啊 ，想学习这方面呢

marlboro_25214

能不能普及一下 SMP MPP什么意思