- 积分
- 0
- 注册时间
- 2011-10-3
- 仿真币
-
- 最后登录
- 1970-1-1
|
本帖最后由 五星连珠 于 2011-11-12 09:17 编辑
开始是在论坛上下的文件。。。。转换的元素有位移、主应力、正应力等 但是我现在需要画一个第二主应力减去第三主应力的量的等值线图。。。我把程序做了一部分改动 但是我的程序总出错 我水平有限 看不出来哪里错了 找人帮解答 最好给我个能达到要求的改后程序,,,,如果你要能给我解释解释就更感谢了
程序如下:
;; by dynamax@simwe
;; 9/14/2004
;; Modifications:
;; (1) change data packing mode from point to block
;; (2) add zone-related info to be plotted
;;
;; 9/11/2004
;; Modifications:
;; (1) add range tec_range to set plot range
;; (2) add function plot_test to hide null zones
;;
;; 7/11/2004
;; Flac3d Mesh to Tecplot
;; Limits:
;; (1)Zone Type ---- TETRAHEDRON or BRICK
;; (2)Tecplot Version ---- Version 10
;; (3)Only DISPLACEMENTS of Grid Points are outputted
;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; Initialization
def ini_mesh2tec
IO_READ = 0
IO_WRITE = 1
IO_FISH = 0
IO_ASCII = 1
N_RECORD = 10
ZONE_NGP = z_numgp(zone_head)
array buf(1)
tec_file = 'tec10.dat'
;; Edit the tec_range to set plot range
command
ran name tec_range
endcommand
end
ini_mesh2tec
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; If plotit = 1, plot the zone
def plot_test
plotit = 1
if z_model(p_z) = 'NULL' then
plotit = 0
endif
if inrange('tec_range',p_z) = 0 then
plotit = 0
endif
end
;; Get number of zones to plot
def get_nzone
n_zone = 0
p_z = zone_head
loop while p_z # null
plot_test
if plotit = 1 then
n_zone = n_zone + 1
endif
p_z = z_next(p_z)
endloop
end
get_nzone
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; Write Tecplot File Head
def write_head
buf(1) = 'TITLE = "FLAC3D to Tecplot 10"\n'
buf(1) = buf(1) + 'VARIABLES = "X(m)" \n"Y(m)" \n"Z(m)" \n'
buf(1) = buf(1) + '"DISP(m)" \n"XDISP(m)" \n"YDISP(m)" \n"ZDISP(m)" \n'
buf(1) = buf(1) + '"SIG1(MPa)" \n"SIG2(MPa)" \n"SIG3(MPa)" \n'
buf(1) = buf(1) + '"SXX(MPa)" \n"SYY(MPa)" \n"SZZ(MPa)" \n'
buf(1) = buf(1) + '"SQ(MPa)"\n'
buf(1) = buf(1) + 'ZONE T="GLOBAL" \n'
buf(1) = buf(1) + ' N=' + string(ngp) + ','
buf(1) = buf(1) + ' E=' + string(n_zone) + ','
if ZONE_NGP = 4 then
buf(1) = buf(1) + ' ZONETYPE=FETETRAHEDRON \n'
else
buf(1) = buf(1) + ' ZONETYPE=FEBrick \n'
endif
buf(1) = buf(1) + ' DATAPACKING=BLOCK \n'
buf(1) = buf(1) + ' VARLOCATION=([8-13]=CELLCENTERED) \n'
buf(1) = buf(1) + ' DT=(SINGLE SINGLE SINGLE'
buf(1) = buf(1) + ' SINGLE SINGLE SINGLE SINGLE'
buf(1) = buf(1) + ' SINGLE SINGLE SINGLE SINGLE'
buf(1) = buf(1) + ' SINGLE SINGLE SINGLE'
buf(1) = buf(1) + ' SINGLE)'
status = write(buf,1)
end
;; Calculate displacement magnitude
def get_gp_disp
gp_disp = gp_xdisp(p_gp)*gp_xdisp(p_gp)
gp_disp = gp_disp + gp_ydisp(p_gp)*gp_ydisp(p_gp)
gp_disp = gp_disp + gp_zdisp(p_gp)*gp_zdisp(p_gp)
gp_disp = sqrt(gp_disp)
end
;; Write gp-related data,such as Coordinates and Displacements
def write_gp_info
p_gp = gp_head
loop while p_gp # null
buf(1) = ''
loop i(1,N_RECORD)
if p_gp # null then
caseof info_flag
case 0
buf(1) = buf(1) + string(gp_xpos(p_gp)) + ' '
case 1
buf(1) = buf(1) + string(gp_ypos(p_gp)) + ' '
case 2
buf(1) = buf(1) + string(gp_zpos(p_gp)) + ' '
case 4
get_gp_disp
buf(1) = buf(1) + string(gp_disp) + ' '
case 8
buf(1) = buf(1) + string(gp_xdisp(p_gp)) + ' '
case 16
buf(1) = buf(1) + string(gp_ydisp(p_gp)) + ' '
case 32
buf(1) = buf(1) + string(gp_zdisp(p_gp)) + ' '
endcase
p_gp = gp_next(p_gp)
endif
endloop
status = write(buf,1)
endloop
end
def get_z_sq
z_sq = z_sig2(p_z) - z_sig3(p_z)
end
;; Write zone-related data,such as Stresses
def write_zone_info
p_z = zone_head
loop while p_z # null
buf(1) = ''
loop i(1,N_RECORD)
if p_z # null then
plot_test
if plotit = 1 then
caseof info_flag
case 0
buf(1) = buf(1) + string(z_sig1(p_z) / (-1000000.0)) + ' '
case 1
buf(1) = buf(1) + string(z_sig2(p_z) / (-1000000.0)) + ' '
case 2
buf(1) = buf(1) + string(z_sig3(p_z) / (-1000000.0)) + ' '
case 4
buf(1) = buf(1) + string(z_sxx(p_z) / (-1000000.0)) + ' '
case 8
buf(1) = buf(1) + string(z_syy(p_z) / (-1000000.0)) + ' '
case 16
buf(1) = buf(1) + string(z_szz(p_z) / (-1000000.0)) + ' '
case 32
get_z_sq
buf(1) = buf(1) + string(z_sq(p_z) / (-1000000.0)) + ' '
endcase
endif
p_z = z_next(p_z)
endif
endloop
status = write(buf,1)
endloop
end
;; Write Zone Connectivity
def write_zone
p_z = zone_head
loop while p_z # null
buf(1) = ' '
if ZONE_NGP = 4 then
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 1))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 2))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 3))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 4))) + ' '
else
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 1))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 2))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 5))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 3))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 4))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 7))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 8))) + ' '
buf(1) = buf(1) + string(gp_id(z_gp(p_z, 6))) + ' '
endif
plot_test
if plotit = 1 then
status = write(buf,1)
endif
p_z = z_next(p_z)
endloop
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;; Main Function
def mesh2tec
status = close
status = open(tec_file,IO_WRITE,IO_ASCII)
if status = 0 then
write_head
info_flag = 0
write_gp_info
info_flag = 1
write_gp_info
info_flag = 2
write_gp_info
info_flag = 4
write_gp_info
info_flag = 8
write_gp_info
info_flag = 16
write_gp_info
info_flag = 32
write_gp_info
info_flag = 0
write_zone_info
info_flag = 1
write_zone_info
info_flag = 2
write_zone_info
info_flag = 4
write_zone_info
info_flag = 8
write_zone_info
info_flag = 16
write_zone_info
info_flag = 32
write_zone_info
write_zone
status = close
ii = out('Successfully Write ' + string(n_zone) + ' Zone Data Into File ' + tec_file)
else
ii=out('Open File Error! Status = ' + string(status))
endif
end
mesh2tec
在FLAC3D里面读入程序 出现如下错误
|
本帖子中包含更多资源
您需要 登录 才可以下载或查看,没有账号?注册
×
|