图1 开挖卸荷后位移分布
图2 雨水只在裂隙中渗流
图3 压力边界条件和流量边界条件示意图
图4 边界条件所对应的高度-含水率-时间图
我在本科期间从未接触过数值计算，但学过Vb程序设计，利用其编过一些小程序完成过作业。而我第一次接触MatDEM是在我研究生入学前，我的导师李小双老师告诉我的，让我利用好暑期学习一下。当时他说:"国产化在任何领域都会是趋势，这个软件以后绝对会成为主流。"但到研究生入学后，我才真正开始学习MatDEM。
一、跟着书本按部就班学习、定一个小目标当时网上的视频资源还没有现在这么丰富，我刚开始学的时候靠得就是书籍。在看了书前面的基本原理后对离散元法有了一定的了解后，我就按部就班地从书上的建模那章开始看，参照书上的实例代码去软件上敲，这让我想起了是以前编制Vb程序写作业的经历：写一段程序去完成一个目标。
这就让我提起了比较浓厚兴趣，研一的课还是比较多的，在利用课余时间在学会如何初始堆积建模、利用高程数据建模、材料赋值并且对相关代码的功能有了一定的理解，再读了点基于MatDEM发表的论文以后我便摩拳擦掌，想要做出点什么东西试试，也就是给自己定一个要做出来的东西，定一个目标。
恰好我当时正在帮一个其他导师带的学姐做一个关于降雨对尾矿库影响的实验，苦力是我当时的身份，实验一直持续到了快放寒假的时候。我才第一次尝试去基于MatDEM完成自己想达到的目标。
二、动手编程去遇见问题、交几个朋友动手去写程序，遇到问题再解决就是推动自己达到自己天花板的最佳途径。如何去考虑降雨和后方的积水入渗呢?如何去考虑其给构筑物带来的影响？在矿山开挖模拟中去考虑蠕变特性？这涉及到多方面的知识，软件本身只是一个工具，想要把这个工具为你所用变得适合你，改进是必须的步骤。
初级的改进方法就是基于它的现有功能进行二次开发。可是在读了相关论文后并且把说明书中的命令都读了后，我并没有找到我力所能及的方法去模拟实验的过程。所以当时本着假期用来练习的想法，做了简化，我降低了自身期望把目标改成了定义一个浸润面、为浸润面以下不同深度的单元颗粒赋予相应的浮力和渗透力。
在这个达成这几个点过程中遇到了很多问题，和很多同学一样我遇到问题一开始也总想着在群里问，但是显然大家都有自己要做的事，没人总会盯着群。此时交几个同样处于苦恼中的朋友一起讨论就十分有必要，有了朋友的陪伴学习MatDEM将会更加欢乐，我就很庆幸我认识了三峡大学的闫金洲同学。
三、如何解决问题、提什么样的问题MatDEM是采用Matlab编译的独立软件，其可以调用Matlab的基本函数进行二次开发。不管是赋值这样的基本操作，还是编写理论公式函数，以及制作软件窗口，都可以通过二次开发完成。
所以我想大多数同学一定会遇到和有关于矩阵运算、处理等相关的命令，就比如：在模拟矿山开挖过程中会有很多冒落颗粒下坠，这时你不想保留它，那么怎么做到既删除又保留不想删除的矿柱呢？这就要群组号进行筛选并且进行“且和或”处理了，又比如想随机生成点裂隙这时就要用到rand()函数命令了。
所以我们在学习中不可避免地需要用Matlab函数对矩阵进行处理，但我又不想系统得学一遍或者好像之前学过点也不会怎么办？举个例子，在MatDEM中的流固耦合模拟中需要通常是需要对邻居矩阵进行处理的，我想要把每个颗粒的含水率（d.mo.SET.aWC）都根据索引号填入邻居矩阵的相应位置该怎么办？最快的方法当然是理解现成案例的写法了，但是如果没有案例，那就需要观察下要达到我们的目标需要用到的原料是什么。
邻居矩阵：行序号为每个颗粒的标号，后面的填好号码都是与他相邻的颗粒的编号。这时手上还有个记录了所有颗粒的含水率的矩阵，其索引号和颗粒号码编号是一致的，再想下目标是填入活动颗粒和边界条件颗粒的含水率，怎么把这两样原料合成想要的东西呢。
再把这个问题描述下：矩阵A是一个n×1矩阵，记录了部分索引号对应的值，矩阵B是n×n矩阵，其中的不同位置记录了A各索引号，请问有方法把B中索引号所对应的A的值填入B的相应位置吗？把这个问题发到群里或者问朋友或者发CSDN上，我想问得这么具体了一定有人能够解答你了，有人会告诉你B其实叫“索引矩阵”。
因此在解决问题和提问这个环节，一定是要基于自己的思考知道的，实在解决不了再详细描述好问题好让别人容易懂有切入点，那么快速得到回复的几率就会高很多了。
四、不断积累知识，有据可依地慢慢改进除了学习编程方面的知识，最核心的还是要补充模拟内容相关的数理化知识。数值计算软件就是把重要的公式编好封装成了软件，有的包装程度高、包装的内容多，想真得搞定数值模拟还是要和那些复杂的公式和数值计算方法打交道。
我的毕业论文是岩样的蠕变模拟，我就是因为才疏学浅没水平在MatDEM里自己开发蠕变接触模型所以毕业论文才去用PFC做的。引入一些假设也许能让模拟出来的效果看上去像回事，但不能和实验、实际工程对应得上是经不起推敲的。
这需要不断得补充关于这方面的知识，看别人的怎么做得，控制方程如何推导的，说不定就能把别人论证过的方法移植过来做到尽量有理有据了，比如之前一直困扰我的入渗问题，当时交流群里有好几个同学同样在尝试的时候，我自己也独立编程过一个案例,结果看着也挺像那么回事。
但是显然，没有如果理论支撑的话只是样子好看，并不能进行实际应用。但在今年读完一篇论文后感觉有戏了，按捺不住激动的立马进行了相应的程序编写，这才有我现在进一步改进案例的可能，以期待能够应用。如果一时半会解决不了问题，那把数值计算作为一项爱好慢慢修炼吧，在以后的工作中也不断取得进步。
总之，作为广大学习者中的中的一个，在其它同学总结的最佳学习方法的基础上，我想补充得就是：
定一个小目标。动手编程主动遇见问题。交几个愿意和你讨论的朋友。学会如何解决问题，把问题描述得尽量详细以提高问题回复率。把解决没完成数值计算问题当成自己的一个爱好，不断积累模拟内容的相关知识进而对曾经的案例进行慢慢改进。五、案例演示下面进行简单的案例演示。
1、降雨强度测试
渗流部分采用VG模型，水力参数保持饱和含水率 ，残余含水率 ,初始含水率 ，渗透系数 ,拟合参数n=2，m=0.5， 0.00335。降雨强度分别为 情况下的渗流情况，降雨时间保持为4h。
雨强为0.00013cm/s时的水分场
不同雨强入深深度对比
从上图可以看出，入渗深度随着降雨强度的增大而增大，其边界条件是会从流量边界自动转换成压力边界的，但由于该计算方法本身的一些缺陷，如果完全饱和会导致计算收敛很慢。
2、降雨入渗与卸荷耦合案例
降雨强度为0.0013cm/s,降雨时间为24h，水力参数与上文参数保持一致。
采动卸荷和卸荷-渗流耦合作用示意图
降雨与不降雨测点竖向位移对比
从上图可以看出，降雨后对岩体会产生影响，简单得将相关参数进行了强度折减，依据宏微观参数直接的关系对微观参数进行强度折减。饱和颗粒细观参数的折减程度为初始情况下参数的80%。可以发现降雨后测点处竖向位移比起未降雨的要大。
六、我的直播公开课10月27日20时，笔者将在仿真秀官网和APP同步在线直播《基于MatDEM》，诚邀大家莅临交流，如有不当，欢迎批评指正，直播支持反复回看。
报名方式：想参与的小伙伴可以私信
作者：周桃涛 仿真秀平台用户 江西理工大学硕士研究生
声明：原创作品，首发仿真秀，部分图片和内容源自网络，如有不当请联系我们，欢迎分享，禁止私自转载，转载请联系我们。
