Python再入手-03

又过了半年，已经完全忘了当时8月份的python工作了。这半年，先是跑合同，然后张罗出海，完了又搞了两次验收，还做了两次汇报，开了无数的会，忙坏了。

现在，得重新张罗电磁这档子事儿了。下面记录下最近的操作。

1 打开编程环境。

都忘了怎么打开环境了，翻一下以前的博客（这是第三篇），终于想起来了。

打开anaconda，然后激活环境e1（总共装了2个，1个e1,1个e2，忘了哪个能用了，随便激活1个）。然后打开spyder。 

下面是spyder界面：

这是4.2.5版本，python用的3.7

2 打开simpeg文件

装这个anaconda主要是为了跑simpeg，里面最感兴趣的就是频率域三维电磁了。下面就逐步分析一下这个文件。

2.1 引入依赖包

这里的代码是如下。首先来盘点一下：

• 从discretize引入treemesh包，实现树形网格；
• 从discretize.utils引入mkvc和refine_tree_xyz，实现向量组装和细化数网格；
• 从SimPEG.utils引入plot2Ddata和surface2ind_topo，实现二维绘图和地形；
• 引入SimPEG的电磁频率域模块，命名为fdem；
• 引入os；
• 引入numpy，命名为np；
• 引入matplotlib，命名为mpl；
• 引入matplotlib中的pyplot，命名为plt
• 从以下选项尝试：先从pymatsolver中引入Pardiso求解器，命名为Solver，如果发生了ImportError，就执行：从SimPEG引入SuperLU求解器，并命名为Solver。
• 规定变量save_file的值为true。

from discretize import TreeMesh
from discretize.utils import mkvc, refine_tree_xyz

from SimPEG.utils import plot2Ddata, surface2ind_topo
from SimPEG import maps
import SimPEG.electromagnetics.frequency_domain as fdem

import os
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt

try:
    from pymatsolver import Pardiso as Solver
except ImportError:
    from SimPEG import SolverLU as Solver

save_file = True

2.2 定义地形

xx, yy = np.meshgrid(np.linspace(-3000, 3000, 101), np.linspace(-3000, 3000, 101))
zz = np.zeros(np.shape(xx))
topo_xyz = np.c_[mkvc(xx), mkvc(yy), mkvc(zz)]  # 定义(N, 3)numpy数组

这段代码写的很紧凑。

• 首先，利用np中的linspace生成向量，然后作为meshgrid的输入变量，生成xx和yy。
• 之后，利用xx的维度，生成0向量作为zz。
• 最后，用mkvc将矩阵变为向量，然后调用c_生成topo_xyz。

2.3 定义观测装置

frequencies = [100, 500, 2500]
N = 9
xtx, ytx, ztx = np.meshgrid(np.linspace(-200, 200, N), np.linspace(-200, 200, N), [40])
source_locations = np.c_[mkvc(xtx), mkvc(ytx), mkvc(ztx)]
ntx = np