1 | initial version |
With some insights from the following posts:
substitute expression instead of formal function symbol
How do I get an ordered list of a symbolic functions arguments?
Call a function with argument list in python (non-sage link)
I've come up with a function, transform_func_x(expr,f,rels), to change the dependent variables of an expression. Below is the function itself and below that how I've used it (all copy and pasted from a sage notebook). Any suggested improvements are welcome.
def transform_func_x(expr,f,rels):
r"""
Transform the dependent variables of a symbolic function
INPUT:
- ``expr`` - symbolic function. Usually a relational expression
with differential operators such as a partial differential equation.
- ``f`` - expression defining a function (e.g. previously defined
as f = function ('u', z, t))
- ``rels`` - list of relational expressions giving variable
transformations. For example to change y(z, t) to y(Z, T_v), rels
could be [Z == z/H, T_v = c_v * t / H^2]. To change y(z, t) to y(Z, t),
rels could be [Z == z / H]
.. note::
1. Transformations will not be performed for any relationship where the
LHS is already one of the variables in f. For, example if
f = u(Z, t) and rels = [Z == z / H, T_v == c_v * t / H^2] then the
"Z == z / H" transformation will be ignored.
2. After applying this function to expr, you will probably want to
update f by applying transform_func_x to f itself
i.e transform_func_x(f,f,rels)
"""
#Limit#all#lines#to#a#maximum#of#79#characters#################################
y = f.operator()
x = f.operands()
nx = f.number_of_operands()
new_arg_list = x[:]
"""used to transform variables within f to functions
of those variables. e.g. u(z, t) --> u(z / H, t)"""
subs_dict1 = {}
"""used to change name of transformed variables within f.
e.g. u(z / H, t) --> u(Z, t)"""
subs_dict2 = {}
"""used to transform variables outside of f.
e.g. z *u(Z, t)--> z * H * u(Z, t)"""
for v in rels:
if v.left() not in x:
"""i.e. exclude transformations where LHS of
rels is already a variable of f"""
for i in range(nx):
# find which variable the transformation applies to
if x[i] in v.right().variables():
xi=i
break
new_arg_list[xi] = v.right()
subs_dict1[v.right()]=v.left()
subs_dict2[x[xi]]=solve(v,x[xi],solution_dict=true)[_sage_const_0 ][x[xi]]
# make dummy function for transforming variables of f
dummy_arg_list = [var("x{0:d}".format(i)) for i in range(nx)]
dummy_y = symbolic_expression(function(str(y),*new_arg_list)).function(*dummy_arg_list)
return (expr.substitute_function(y,dummy_y).subs(subs_dict1).
subs(subs_dict2))
reset
var('z,t,H,T_v,Z,c_v,k_v,m_v, g_w')
f=function('u',z,t)
pde = k_v*f.diff(z,z)/g_w == m_v*f.diff(t)
print('Variable transformations')
print('PDE with u(z,t):')
rels=[]
eq1=transform_func_x(pde,f,rels)
f=transform_func_x(f,f,rels)
print(eq1)
print('Term to isolate:')
myterm = f.diff(f.operands()[0],2)
print (myterm)
print('Coefficient of term on LHS:')
mycoeff=eq1.left_hand_side().coeff(myterm)
print(mycoeff)
print ('Divide both sides by coefficient:')
eq2 = eq1.divide_both_sides(mycoeff)
print(eq2)
print('Term to isolate:')
myterm = f.diff(f.operands()[1])
print(myterm)
print('Coefficient of term on RHS:')
mycoeff = eq2.right_hand_side().coeff(myterm)
print(mycoeff)
print ('Set RHS coefficient to c_v:')
eq3 = c_v == 1/mycoeff
print(eq3)
print('Rearrange c_v equation for substitution:')
eq4 = solve(eq3,eq3.right_hand_side().default_variable(), solution_dict=true)
print(eq4)
print('Make c_v substitution:')
eq5 = eq2.subs(eq4[0])
print(eq5)
print('\nz is in the domain [0,H] so introduce variable change Z = z/H:')
rels=[Z==z/H]
print(rels)
print('PDE with u(Z,t):')
eq6=transform_func_x(eq5,f,rels)
f=transform_func_x(f,f,rels)
print (eq6)
print('Term to isolate:')
myterm = f.diff(f.operands()[0],2)
print (myterm)
print('Coefficient of term on LHS:')
mycoeff=eq6.left_hand_side().coeff(myterm)
print(mycoeff)
print ('Divide both sides by coefficient:')
eq7 = eq6.divide_both_sides(mycoeff)
print(eq7)
print('\nIntroduce variable change T_v = c_v * t / H^2 to simplify RHS:')
rels=[Z == z/H,T_v==c_v*t/H^2]
print(rels)
print('PDE with u(Z,T_v):')
eq8=transform_func_x(eq7,f,rels)
f=transform_func_x(f,f,rels)
print (eq8)
Running the above I get:
Variable transformations
PDE with u(z,t):
k_v*D[0, 0](u)(z, t)/g_w == m_v*D[1](u)(z, t)
Term to isolate:
D[0, 0](u)(z, t)
Coefficient of term on LHS:
k_v/g_w
Divide both sides by coefficient:
D[0, 0](u)(z, t) == g_w*m_v*D[1](u)(z, t)/k_v
Term to isolate:
D[1](u)(z, t)
Coefficient of term on RHS:
g_w*m_v/k_v
Set RHS coefficient to c_v:
c_v == k_v/(g_w*m_v)
Rearrange c_v equation for substitution:
[{g_w: k_v/(c_v*m_v)}]
Make c_v substitution:
D[0, 0](u)(z, t) == D[1](u)(z, t)/c_v
z is in the domain [0,H] so introduce variable change Z = z/H:
[Z == z/H]
PDE with u(Z,t):
D[0, 0](u)(Z, t)/H^2 == D[1](u)(Z, t)/c_v
Term to isolate:
D[0, 0](u)(Z, t)
Coefficient of term on LHS:
H^(-2)
Divide both sides by coefficient:
D[0, 0](u)(Z, t) == H^2*D[1](u)(Z, t)/c_v
Introduce variable change T_v = c_v * t / H^2 to simplify RHS:
[Z == z/H, T_v == c_v*t/H^2]
PDE with u(Z,T_v):
D[0, 0](u)(Z, T_v) == D[1](u)(Z, T_v)
2 | No.2 Revision |
With some insights from the following posts:
substitute expression instead of formal function symbol
How do I get an ordered list of a symbolic functions arguments?
Call a function with argument list in python (non-sage link)
I've come up with a function, transform_func_x(expr,f,rels), to change the dependent variables of an expression. Below is the function itself and below that how I've used it (all copy and pasted from a sage notebook). Any suggested improvements are welcome.
def transform_func_x(expr,f,rels):
r"""
Transform the dependent variables of a symbolic function
INPUT:
- ``expr`` - symbolic function. Usually a relational expression
with differential operators such as a partial differential equation.
- ``f`` - expression defining a function (e.g. previously defined
as f = function ('u', z, t))
- ``rels`` - list of relational expressions giving variable
transformations. For example to change y(z, t) to y(Z, T_v), rels
could be [Z == z/H, T_v = c_v * t / H^2]. To change y(z, t) to y(Z, t),
rels could be [Z == z / H]
.. note::
1. Transformations will not be performed for any relationship where the
LHS is already one of the variables in f. For, example if
f = u(Z, t) and rels = [Z == z / H, T_v == c_v * t / H^2] then the
"Z == z / H" transformation will be ignored.
2. After applying this function to expr, you will probably want to
update f by applying transform_func_x to f itself
i.e transform_func_x(f,f,rels)
"""
#Limit#all#lines#to#a#maximum#of#79#characters#################################
y = f.operator()
x = f.operands()
nx = f.number_of_operands()
new_arg_list = x[:]
"""used to transform variables within f to functions
of those variables. e.g. u(z, t) --> u(z / H, t)"""
subs_dict1 = {}
"""used to change name of transformed variables within f.
e.g. u(z / H, t) --> u(Z, t)"""
subs_dict2 = {}
"""used to transform variables outside of f.
e.g. z *u(Z, t)--> z * H * u(Z, t)"""
for v in rels:
if v.left() not in x:
"""i.e. exclude transformations where LHS of
rels is already a variable of f"""
for i in range(nx):
# find which variable the transformation applies to
if x[i] in v.right().variables():
xi=i
break
new_arg_list[xi] = v.right()
subs_dict1[v.right()]=v.left()
subs_dict2[x[xi]]=solve(v,x[xi],solution_dict=true)[_sage_const_0 ][x[xi]]
# make dummy function for transforming variables of f
dummy_arg_list = [var("x{0:d}".format(i)) for i in range(nx)]
dummy_y = symbolic_expression(function(str(y),*new_arg_list)).function(*dummy_arg_list)
dummy_y(*dummy_arg_list) = function(str(y),*new_arg_list))
return (expr.substitute_function(y,dummy_y).subs(subs_dict1).
subs(subs_dict2))
reset
var('z,t,H,T_v,Z,c_v,k_v,m_v, g_w')
f=function('u',z,t)
pde = k_v*f.diff(z,z)/g_w == m_v*f.diff(t)
print('Variable transformations')
print('PDE with u(z,t):')
rels=[]
eq1=transform_func_x(pde,f,rels)
f=transform_func_x(f,f,rels)
print(eq1)
print('Term to isolate:')
myterm = f.diff(f.operands()[0],2)
print (myterm)
print('Coefficient of term on LHS:')
mycoeff=eq1.left_hand_side().coeff(myterm)
print(mycoeff)
print ('Divide both sides by coefficient:')
eq2 = eq1.divide_both_sides(mycoeff)
print(eq2)
print('Term to isolate:')
myterm = f.diff(f.operands()[1])
print(myterm)
print('Coefficient of term on RHS:')
mycoeff = eq2.right_hand_side().coeff(myterm)
print(mycoeff)
print ('Set RHS coefficient to c_v:')
eq3 = c_v == 1/mycoeff
print(eq3)
print('Rearrange c_v equation for substitution:')
eq4 = solve(eq3,eq3.right_hand_side().default_variable(), solution_dict=true)
print(eq4)
print('Make c_v substitution:')
eq5 = eq2.subs(eq4[0])
print(eq5)
print('\nz is in the domain [0,H] so introduce variable change Z = z/H:')
rels=[Z==z/H]
print(rels)
print('PDE with u(Z,t):')
eq6=transform_func_x(eq5,f,rels)
f=transform_func_x(f,f,rels)
print (eq6)
print('Term to isolate:')
myterm = f.diff(f.operands()[0],2)
print (myterm)
print('Coefficient of term on LHS:')
mycoeff=eq6.left_hand_side().coeff(myterm)
print(mycoeff)
print ('Divide both sides by coefficient:')
eq7 = eq6.divide_both_sides(mycoeff)
print(eq7)
print('\nIntroduce variable change T_v = c_v * t / H^2 to simplify RHS:')
rels=[Z == z/H,T_v==c_v*t/H^2]
print(rels)
print('PDE with u(Z,T_v):')
eq8=transform_func_x(eq7,f,rels)
f=transform_func_x(f,f,rels)
print (eq8)
Running the above I get:
Variable transformations
PDE with u(z,t):
k_v*D[0, 0](u)(z, t)/g_w == m_v*D[1](u)(z, t)
Term to isolate:
D[0, 0](u)(z, t)
Coefficient of term on LHS:
k_v/g_w
Divide both sides by coefficient:
D[0, 0](u)(z, t) == g_w*m_v*D[1](u)(z, t)/k_v
Term to isolate:
D[1](u)(z, t)
Coefficient of term on RHS:
g_w*m_v/k_v
Set RHS coefficient to c_v:
c_v == k_v/(g_w*m_v)
Rearrange c_v equation for substitution:
[{g_w: k_v/(c_v*m_v)}]
Make c_v substitution:
D[0, 0](u)(z, t) == D[1](u)(z, t)/c_v
z is in the domain [0,H] so introduce variable change Z = z/H:
[Z == z/H]
PDE with u(Z,t):
D[0, 0](u)(Z, t)/H^2 == D[1](u)(Z, t)/c_v
Term to isolate:
D[0, 0](u)(Z, t)
Coefficient of term on LHS:
H^(-2)
Divide both sides by coefficient:
D[0, 0](u)(Z, t) == H^2*D[1](u)(Z, t)/c_v
Introduce variable change T_v = c_v * t / H^2 to simplify RHS:
[Z == z/H, T_v == c_v*t/H^2]
PDE with u(Z,T_v):
D[0, 0](u)(Z, T_v) == D[1](u)(Z, T_v)