首页/文章/ 详情

基于Python脚本提取黄永刚晶体塑性miller指数并转化为欧拉角绘制极图

2月前浏览2307

前述介绍了几种方案用于提取黄永刚晶体塑性的欧拉角文件绘制极图,如

(1)在huang子程序中进行晶体取向更新 (qq.com)------即修改umat问价加入状态变量提取

(2)基于lingzhi-matlab脚本与huang-umat实现变形过程中的织构演化预测-技 术邻 (jishulink.com) ------在for中写出miller指数并使用matlab转化为欧拉角绘制极图

这两类方案均可以较为方便的绘制极图,但方案一个别情况下会导致出现计算中的矩阵奇异,而方式二计算过程中因为线程原因无法并行运算,此外两类方案先显示中实现的效果通常极不方便。

因此这里介绍第三类方法绘制黄永刚后处理的极图-------使用Python脚本进行状态变量的提取,然后使用matlab直接绘制极图,这类方案的的主要优势是无需修改对应的Fortran程序即可实现,对原始程序的收敛性无影响,而且适用于显示vumat

程序来源于技术邻icpfem的开源分享:晶体塑性有限元仿真入门(5)—欧拉角与晶体取向_MATLAB ABAQUS-技 术邻 (jishulink.com)

简单介绍如下:

(1)使用黄永刚程序计算完之后,利用Python脚本提取miller指数:主要脚本如下(展示部分内容,完整代码可以查看原始网页):

# SDV提取,本程序适用于将odb输出文件的SDV提取导出保存

     

import odbAccess

import numpy as np

     

myOdb = odbAccess.openOdb(r"D:\CPFEM_Abaqus_V2\Course3_3\Job-1201-001.odb",readOnly=False)

myStep = myOdb.steps["Step-1"]

myFrame = myStep.frames

i = len(myFrame)-1

myField = myFrame[i].fieldOutputs

AA = 3;

myValue = myField["SDV37"].values

temp_h = []

for value in myValue:  

   myData = value.data

   temp_h.append(myData)

  ..........   


(2)之后将miller指数数据转化为euler角并存储,这一步使用matlab完成

clc

clear

close all

% Data = load('D:\CPFEM_Abaqus_V2\Course3_3\Out_end_rolling.txt');

Data = load('D:\CPFEM_Abaqus_V2\Course3_3\Out_end_compression.txt');

HKLandUVW = reshape(Data,length(Data)/6,6);

data = HKLandUVW(1:1:end,:);

Unique_data = data(1:20:end,:);

Euler = zeros(length(Unique_data(:,1)),3);

for i = 1:length(Unique_data(:,1))

   Euler(i,:) = SDV_to_Euler(Unique_data(i,1:3),Unique_data(i,4:6));

............

(3)利用转化得到的欧拉角绘制极图

cs = crystalSymmetry('cubic');        

ss = specimenSymmetry('triclinic');

ori = orientation.byEuler(Euler(:,1)*degree,Euler(:,2)*degree,Euler(:,3)*degree,cs,ss);

plotPDF(ori,Miller({1,0,0},{1,1,0},{1,1,1},cs),'all','MarkerSize',1)

pause(5)

     

Index = transpose(1:length(Euler(:,1)));

Phase = ones(length(Euler(:,1)),1);

X = transpose(1:length(Euler(:,1)));

Y = transpose(1:length(Euler(:,1)));

MAD = rand(length(Euler(:,1)),1);

BC = ceil(100*rand(length(Euler(:,1)),1));

BS = zeros(length(Euler(:,1)),1);

..............

这里展示作者验证的结果

可以看到使用这种方法可以很方便且正确的绘制极图,同时我使用该方案如前述方案进行了比较,发现结果具有良好的一致性,但原始网页未进行显示vumat程序的对比

这里对比使用前述方案与Python提取方法在显示分析中晶体织构绘制的差异,验证模型在显示的适用性。

多晶拉伸结束后的极图(显示结果)

(a)以前的方案与python提取状态变量方案绘制极图分别如下

可以看到两者具有很好的一致性。因此显示分析中该脚本同样适用。感兴趣的可以查看原始网页查看相关脚本。使用脚本过程中存在问题可以加入我的知识星球进行讨论交流,另外脚本使用的完整流程也可以加入我的知识星球进行获取。




来源:我的博士日记
AbaqusSTEPSMATLABpythonUM
著作权归作者所有,欢迎分享,未经许可,不得转载
首次发布时间:2024-09-10
最近编辑:2月前
此生君子意逍遥
博士 签名征集中
获赞 48粉丝 66文章 84课程 0
点赞
收藏
作者推荐

超弹性与亚弹性,显式与隐式的HCP多晶滑移+孪晶(主导孪晶重定向(PTR))计算效率比较

参考文献:《Influence of texture distribution in magnesium welds on their non-uniform mechanical behavior: A CPFEM study》主导孪晶重定向(PTR)方案作为目前处理HCP晶格结构的多晶材料孪晶模拟中最常使用的方案被广泛讨论,然而晶体取向旋转过程可能会造成模拟的收敛性问题,选择一个相对稳定的本构框以及迭代变量对模拟计算效率的提升是有意义的。文章中公式以及其对应的参数总结如下:这里使用文章的模型和参数对超弹性和亚弹性PTR方案进行比较。二维(200个晶粒X方向压缩20%)以下各个图中左图为超弹性结果,右图为亚弹性结果:应力分布云图应变分布云图:孪晶分布云图:这里使用文章的模型和参数对显示和隐式PTR方案进行比较二维(200个晶粒Y方向剪切变形20%)以下各个图中左图为显示结果,右图为隐式结果:应力分布云图应变分布云图:孪晶分布云图拉压非对称与织构演化方面超弹性与亚弹性保持一致:初始极图:RD拉伸20%:RD压缩20%:应力应变曲线模拟的结果建议,使用PTR方案,超弹性建议使用PK2应力和当前强度为迭代变量,并使用双重迭代方案,亚弹性建议使用柯西应力为迭代变量,两者在模拟过程中,计算效率相差较小,无论是局部晶粒的应力应变响应,整体的流动应力,以及变形后的织构结果几乎保持一致。同时涉及到接触,碰撞问题,修改为显式对于收敛性的提升是必要的。来源:我的博士日记

未登录
1条评论
仿真秀1014205933
签名征集中
1月前
您好,我在尝试使用该程序以后有几个问题想请教一下,该程序是不是只能绘制整个模型的取向分布结果呢?绘制不同的增量步下的织构情况是通过控制myframe[i]中的i值吗?如果我想绘制某一个截面的取向分布可以吗?
回复
课程
培训
服务
行家
VIP会员 学习 福利任务 兑换礼品
下载APP
联系我们
帮助与反馈