Sfepyのインストール

追記：よくわからんけど色々やっていたら動いた。
最終的には export ETS_TOOLKIT="qt" が鍵だった。と思ったら今度はこの環境変数があると動かなくなった(unsetして消した)。謎が多い。
とりあえず色々エラーの原因になるようなのでAnaconda環境をすべて削除した。

Python 3.8.5 (default, Sep  4 2020, 02:22:02)
[Clang 10.0.0 ] :: Anaconda, Inc. on darwin
Type "help", "copyright", "credits" or "license" for more information.
>>> from mayavi import mlab
multiple 'pyface.toolkits' plugins found for toolkit 'qt': pyface.ui.qt4.init, pyface.ui.qt4.init
multiple 'tvtk.toolkits' plugins found for toolkit 'qt': tvtk.pyface.ui.qt4.init, tvtk.pyface.ui.qt4.init
>>> mlab.plot_3d()
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
AttributeError: module 'mayavi.mlab' has no attribute 'plot_3d'
>>> mlab.test_plot3d()
qt.qpa.window: <QNSWindow: 0x7fb661cfac00; contentView=<QNSView: 0x7fb661cfa3a0; QCocoaWindow(0x7fb661cfa290, window=QWidgetWindow(0x7fb661cf9de0, name="QMainWindowClassWindow"))>> has active key-value observers (KVO)! These will stop working now that the window is recreated, and will result in exceptions when the observers are removed. Break in QCocoaWindow::recreateWindowIfNeeded to debug.
<mayavi.modules.surface.Surface object at 0x7fb6d9e48180>

github.com

Sfepyの可視化に必要なMayaviのインストールが鬼門だった。
Mayavi v4.7.3はvtk9.0.1にバージョンを合わせる必要がある。
PyQt5は5.12.1を指定する。

docs.enthought.com
github.com

$ conda update conda anaconda
$ conda install -c conda-forge sfepy
$ pip install vtk=9.0.1
$ pip install PyQt5==5.12.1 pyqtwebengine==5.12.1
$ pip install mayavi=4.7.3 -vvv

他の依存パッケージ

cython
numpy
scipy
matplotlib
pyparsing
tables
sympy
igakit
petsc4py
mpi4py
slepc4py
pymetis
scikits.umfpack
meshio
psutil
pyvista
opt_einsum
jax
dask

Sfepyのマニュアルにサンプルコードがある。
https://sfepy.org/doc-devel/_downloads/63e3f3e0a0635ddc67af3c6e4d0f6b5f/sfepy_manual.pdf

import numpy as nm
from sfepy.discrete.fem import Mesh, FEDomain, Field
from sfepy import data_dir
filename_mesh = data_dir + '/meshes/2d/rectangle_tri.mesh'
mesh = Mesh.from_file(filename_mesh)
domain = FEDomain('domain', mesh)
min_x, max_x = domain.get_mesh_bounding_box()[:, 0]
eps = 1e-8 * (max_x - min_x)
omega = domain.create_region('Omega', 'all')
gamma1 = domain.create_region('Gamma1','vertices in x < %.10f' % (min_x + eps),'facet')
gamma2 = domain.create_region('Gamma2','vertices in x > %.10f' % (max_x - eps),'facet')                              
field = Field.from_args('fu', nm.float64, 'vector', omega, approx_order=2)
from sfepy.discrete import (FieldVariable, Material, Integral, Function, Equation, Equations, Problem)
u = FieldVariable('u', 'unknown', field)
v = FieldVariable('v', 'test', field, primary_var_name='u')
from sfepy.mechanics.matcoefs import stiffness_from_lame
m = Material('m', D=stiffness_from_lame(dim=2, lam=1.0, mu=1.0)) 
f = Material('f', val=[[0.02], [0.01]])
integral = Integral('i', order=3)
from sfepy.terms import Term
t1 = Term.new('dw_lin_elastic(m.D, v, u)',integral, omega, m=m, v=v, u=u)
t2 = Term.new('dw_volume_lvf(f.val, v)', integral, omega, f=f, v=v)
eq = Equation('balance', t1 + t2)
eqs = Equations([eq])
from sfepy.discrete.conditions import Conditions, EssentialBC
fix_u = EssentialBC('fix_u', gamma1, {'u.all' : 0.0})
def shift_u_fun(ts, coors, bc=None, problem=None, shift=0.0):
    val = shift * coors[:,1]**2
    return val
bc_fun = Function('shift_u_fun', shift_u_fun,extra_args={'shift' : 0.01})
shift_u = EssentialBC('shift_u', gamma2, {'u.0' : bc_fun})
from sfepy.base.base import IndexedStruct 
from sfepy.solvers.ls import ScipyDirect 
from sfepy.solvers.nls import Newton
ls = ScipyDirect({})
nls_status = IndexedStruct()
nls = Newton({}, lin_solver=ls, status=nls_status)
pb = Problem('elasticity', equations=eqs)
pb.save_regions_as_groups('regions')
from sfepy.postprocess.viewer import Viewer
view = Viewer('regions.vtk')
view()

実行結果