CAD转GOCAD转SUFER转FLAC3D建立的复杂地质体模型
我参考freddie的帖子，建立的模型，原帖见http://www.simwe.com/forum/viewthread.php?tid=528370&highlight=%CE%D2%D2%B2%C0%B4surfer%2Bto%2Bflac命令流如下
请教大家：这个是一个层面的，均质体，我对FISH的学习太少，哪位能帮我改下，建立个多层面的模型，谢谢
我的具体方法过程见我的另外帖子http://www.simwe.com/forum/thread-829856-1-5.html
我把命令流贴下，大家帮我看看如何的能改写FISH 建立2层以及多层的地质体模型
主程序
;===============
;main.txt
;===============
new
;====================================================================================
;首先确定网格参数
def pre
oo=out('===============Please input the parameter of the mesh===============')
oo=out('=============== attention ===============')
oo=out('=============== Index_x:Gridpoints along the X-AXE ===============')
Index_x=in('Index_x=')
oo=out('=============== Index_y:Gridpoints along the Y-AXE ===============')
Index_y=in('Index_y=')
oo=out('=============== Index_z:Gridpoints along the Z-AXE ===============')
Index_z=in('Index_z=')
oo=out('=============== Hb_surf:z_cordinate of the below surf ==============')
Hb_surf=in('Hb_surf=')
end
pre
;============= note ============
;实际沿x轴生成的单元数为Index_x-1
;实际沿y轴生成的单元数为Index_y-1
;实际沿z轴生成的单元数为Index_z-1
;===================================
;====================================================================================
;====================================================================================
;调用fish段SURFER.fis,该段得到surfer程序生成的面网格
;用到的参数:Index_x,Index_y
;返回的参数:数组surf_nd存储面网格节点坐标，类型为字符串,状态变量key,节点循环变量loopid
;维数：(Index_x×Index_y)×1×1
call SURFER.fis
;====================================================================================
;====================================================================================
;调用fish段NODEGEN.fis,该段得到surfer程序生成的面网格
;用到的参数:Index_x,Index_y,Index_z,Hb_surf,surf_nd,key
;返回的参数:三维数组node_cor存储体网格节点坐标
;维数：(Index_x×Index_y)×3×Index_z
call NODEGEN.fis
;====================================================================================
;====================================================================================
;调用fish段MESHGEN.fis,该段得到体网格
;用到的参数:Index_x,Index_y,Index_z,node_cor,key
call MESHGEN.fis
;====================================================================================
save mesh.sav
;==========================================================
;SURFER.fis
;用到的参数:Index_x,Index_y
;返回的参数:数组surf_nd存储面网格节点坐标，类型为字符串
;维数：(Index_x×Index_y)×1×1
;==========================================================
;==========================================================
;初始化
def setup
nd_size=Index_x*Index_y ;面节点的总数，也是数组surf_nd的维数
Io_read=0 ;文件形式:surfnd.dat必须存在
I0_ascll=1 ;以ASCII码形式读取文件
filename='surfnd.dat' ;surfer程序产生的面节点文件名
end
setup
;==========================================================
def snd_read
array surf_nd (nd_size)
f_io=open(filename,Io_read,I0_ascll)
status=read(surf_nd,nd_size)
caseof abs(status)
case 1
oo=out('Fatal error:Unexcept end-of-file,you shoult rebuild the surfnd.dat!')
key=1
exit
case 0
00=out('Reading is well done!')
loopid=nd_size
endcase
if abs(status) > 1
oo=out('Fatal error:surf_nd demensions are too large!')
key=1
exit
endif
f_close=close
end
snd_read
;====================================================================================
;调用fish段NODEGEN.fis,该段得到surfer程序生成的面网格
;用到的参数:Index_x,Index_y,Index_z,Hb_surf,surf_nd
;返回的参数:三维数组node_cor存储体网格节点坐标
;维数：(Index_x×Index_y)×3×Index_z
def ndgen
if key=1
exit
endif
arraynode_cor (loopid,3,Index_z)
;z_depth=float(Hb_surf)/(Index_z-1) ;z方向
loop i (1,Index_z)
loop ii (1,loopid)
node_cor(ii,1,i)=parse(surf_nd(ii), 1)
node_cor(ii,2,i)=parse(surf_nd(ii), 2)
Height=(Index_z-i)*(parse(surf_nd(ii), 3)-float(Hb_surf))/(Index_z-1)
node_cor(ii,3,i)=parse(surf_nd(ii), 3)-Height
endloop
endloop
end
ndgen
;====================================================================================
;调用fish段MESHGEN.fis,该段得到体网格
;用到的参数:Index_x,Index_y,Index_z,node_cor,key
;====================================================================================
def mshgen
if key=1
exit
endif
array p1loc (3) p2loc (3) p3loc (3) p4loc (3)
array p5loc (3) p6loc (3) p7loc (3) p8loc (3)
array gpid(4)
loop i (1,Index_z-1)
loop ii (1,Index_x-1)
loop iii (1,Index_y-1)
gpid(1)=(iii-1)*Index_x+ii
gpid(2)=gpid(1)+1
gpid(3)=gpid(1)+Index_x+1
gpid(4)=gpid(3)-1
loop mshid (1,3)
p1loc(mshid)=node_cor(gpid(1),mshid,i)
p2loc(mshid)=node_cor(gpid(2),mshid,i)
p3loc(mshid)=node_cor(gpid(3),mshid,i)
p4loc(mshid)=node_cor(gpid(4),mshid,i)
p5loc(mshid)=node_cor(gpid(1),mshid,i+1)
p6loc(mshid)=node_cor(gpid(2),mshid,i+1)
p7loc(mshid)=node_cor(gpid(3),mshid,i+1)
p8loc(mshid)=node_cor(gpid(4),mshid,i+1)
endloop
p0x=p1loc(1)p0x p0y p0z
p0y=p1loc(2)
p0z=p1loc(3)
p1x=p2loc(1)
p1y=p2loc(2)
p1z=p2loc(3)
p2x=p4loc(1)
p2y=p4loc(2)
p2z=p4loc(3)
p3x=p5loc(1)
p3y=p5loc(2)
p3z=p5loc(3)
p4x=p3loc(1)
p4y=p3loc(2)
p4z=p3loc(3)
p5x=p8loc(1)
p5y=p8loc(2)
p5z=p8loc(3)
p6x=p6loc(1)
p6y=p6loc(2)
p6z=p6loc(3)
p7x=p7loc(1)
p7y=p7loc(2)
p7z=p7loc(3)
command
gen zone brick size 1 1 1 p0 p0x p0y p0z &
p1 p1x p1y p1z p2 p2x p2y p2z) &
p3 p3x p3y p3z p4 p4x p4y p4z &
p5 p5x p5y p5z p6 p6x p6y p6z &
p7 p7x p7y p7z
endcommand
endloop
endloop
endloop
end
mshgen
[ 本帖最后由 青羽儿 于 2008-5-11 20:52 编辑 ]
回复 9# 的帖子
请赐教，是哪出了问题 CAD图[ 本帖最后由 青羽儿 于 2008-5-15 15:23 编辑 ] gocad的地质图
[ 本帖最后由 青羽儿 于 2008-5-15 15:23 编辑 ] flac3d的地质图
[ 本帖最后由 青羽儿 于 2008-5-15 15:23 编辑 ] 请教 ：我用海棠的那个帖子建立的模型地表起伏效果并不好，但是他的FISH可以划分多层，请教如何能结合二者建立多层体 均质体用TOP.FIS就可以建立单层，地面起伏效果是不错的。 斑竹加油啊
回复 6# 的帖子
大哥您竟然我回我的帖子谢谢请教您两个问题1，2层的我用海棠的老是出问题，地表起伏不好，我想用freddie的程序改，可是我没思路，您提示下好吧
2，看了您的那个精华帖子，关于水面的处理（我刚弄这个，还没考虑水面问题），等我弄到这一布的时候就参考您的是不是就行，还用做什么修改么，另外如果变化的水面要怎末考虑呢
谢谢仰慕您的大名 个人人为海棠大哥的程序在形成flac模型时循环少一次，造成边缘部分生成网格的时候会发生错误 海棠大哥那个确实少循环了一次啊！
网格点和网格数目就差了个1！:lol
回复10# 的帖子
海棠兄的03_Gen_FlacModel.fis文件中的zr1=table(j+n_zon_row,x1)
zr2=table(j+n_zon_row+1,x1)
zr3=table(j+n_zon_row+1,x2)
zr4=table(j+n_zon_row,x2)
似乎要改成：
zr1=table(j+Row_Num,x1)
zr2=table(j+Row_Num+1,x1)
zr3=table(j+Row_Num+1,x2)
zr4=table(j+Row_Num,x2) 恩！对就是这样的，开始没有发现是因为海棠大哥的数据刚好能运行！
后来数据换了就出现错误，才发现这个错误！
地震了！慌啊！
回复 12# 的帖子
谢谢哈 哥们 问你个问题 、：如何给模型加上水面回复 9# jianjun_hu 的帖子
是啊，冷风说的对，我也发现有这个问题！！！ 试验了一下修改的，好像出现了问题！ 两层之间出现了大的裂缝，有些地方网格严重变形，局部还多出了网格
页:
[1]