验证ModelicasmoothOrder注释的正确性
这是我之前的问题的后续问题:
三次多项式的 Modelica smoothOrder
我正在为一个大型现有库做出贡献,函数中可能有 50 个 smoothOrder 注释。这些通常是分段定义的函数,在 if-else 语句的不同分支中具有不同的函数定义。
检查平滑度的手动方法可以是评估和绘制接近 if-else 切换条件的函数和所有(偏)导数,直至定义的阶数。如果导数的值和图是连续的,则注释是正确的。
但是这种手动方法非常耗时,因此最好有某种自动检查。这是否存在,或者任何人都可以共享示例脚本或模型来帮助入门?
回答
我已经使用 sympy 或多或少地做了你所描述的。
下面的 jupyter 笔记本计算所有输入方程关于给定变量的偏导数,直到给定的阶数。然后它评估特定值的导数(if 条件中的切换点)。如果输出值相等,则该导数的方程应该是平滑的。
# In[1]:
from IPython.display import display
from sympy import init_printing, symbols, diff, simplify
init_printing()
# In[2]:
def der_at(eqs, wrt, nder, at):
""" Display derivatives of equations up to given order and evaluate for given value """
ordinal = lambda n: "%d%s" % (n,"tsnrhtdd"[(n//10%10!=1)*(n%10<4)*n%10::4])
for i in range(0, nder+1):
der_eqs = [diff(eq, wrt, i) for eq in eqs]
der_eqs_subs = [simplify(der_eq.subs([(wrt, at)])) for der_eq in der_eqs]
print(f"{ordinal(i)}-order derivatives:")
display(der_eqs)
print(f"Derivatives evaluated for {wrt}={at}")
display(der_eqs_subs)
# In[3]:
a, b, c, d, x = symbols('a b c d x')
eq1 = a*x**2 + b*x
eq2 = c*x**3 + d*x**2
der_at(eqs=(eq1, eq2), wrt=x, nder=3, at=0)
运行上面的代码后,您应该得到图像中显示的输出。对于我使用的示例方程,您可以看到,当我们评估 的原始方程时,我们得到了相同的结果x=0,但不是一阶导数。因此平滑阶数为 0。