Solver Issue No Algorithm Found

Viewed 61

Trying to replace the function U2(x,y,z) with specified values of x,y,z. Not sure how to do that with sympy because they are as "x = arange.(-h,h,0.001)" as seen in the code below. Below you will find my implementation with sympy. Additionally I am using PyCharm.

This implementation is based on Dr. Annabestani and Dr. Naghavis' paper: A 3D analytical ion transport model for ionic polymer metal composite actuators in large bending deformations

import sympy as sp
h = 0.1  # [mm] half of thickness
W: float = 6  # [mm] width
L: float = 28  # [mm] length
F: float = 96458  # [C/mol] Faraday's constant
k_e = 1.34E-6  # [F/m]
Y = 5.71E8  # [Pa]
d = 1.03 - 11  # [m^2/s] diffiusitivity coefficient
T = 293  # [K]
C_minus = 1200  # [mol/m^3] Cation concentration
C_plus = 1200 # [mol/m^3] anion concentration
R = 8.3143  # [J/mol*K] Gas constant
Vol = 2*h*W*L
#dVol = diff(Vol,x) + diff(Vol, y) + diff(Vol, z)  # change in Volume
theta = 1 / W
x, y, z, m, n, p, t = sp.symbols('x y z m n p t')
V_1 = 0.5 * sp.sin(2 * sp.pi * t)  # Voltage as a function of time
k_f = 0.5
t_f = 44
k_g = 4.5
t_g = 0.07
B_mnp = 0.003
b_mnp: float = B_mnp
gamma_hat_2 = 0.04
gamma_hat_5 = 0.03
gamma_hat_6 = 5E-3
r_M = 0.15  # membrane resistance
r_ew = 0.175  # transverse resistance of electrode
r_el = 0.11  # longitudinal resistance of electrode
mu = 2.4
sigma_not = 0.1
a_L: float = 1.0  # distrubuted surface attentuation
r_hat = sp.sqrt(r_M ** 2 + r_ew ** 2 + r_el ** 2)
lambda_1 = 0.0001
dVol = 1
K = (F ** 2 * C_minus * d * (1 - C_minus * dVol)) / (R * T * k_e)  # also K = a
K_hat = (K-lambda_1)/d
gamma_1 = 1.0
gamma_2 = 1.0
gamma_3 = 1.0
gamma_4 = 1.0
gamma_5 = 1.0
gamma_6 = 1.0
gamma_7 = 1.0
small_gamma_1 = 1.0
small_gamma_2 = 1.0
small_gamma_3 = 1.0

psi = gamma_1*x + gamma_2*y + gamma_3*z + gamma_4*x*y + gamma_5*x*z + gamma_6*y*z + gamma_7*x*y*z + (small_gamma_1/2)*x**2 + (small_gamma_2/2)*y**2 + (small_gamma_3/2)*x*z**2

psi_hat_part = ((sp.sin(((m + 1) * sp.pi) / 2 * h)) * x) * ((sp.sin(((n + 1) * sp.pi) / W)) * y) * ((sp.sin(((p + 1) * sp.pi) / L)) * z)

psi_hat = psi * psi_hat_part  # Eqn. 19
print(psi_hat)

x1: float = -h
x2: float = h
y1: float = 0
y2: float = W
z1: float = 0
z2: float = L

I = psi_hat.integrate((x, x1, x2), (y, y1, y2), (z, z1, z2))  # Integration for a_mnp Eqn. 18

A_mnp = ((8 * K_hat) / (2 * h * W * L)) * I
Partial = A_mnp * ((sp.sin(((m + 1) * sp.pi) / 2 * h)) * x) * ((sp.sin(((n + 1) * sp.pi) / W)) * y) * ((sp.sin(((p + 1) * sp.pi) / L)) * z)
start = Partial.integrate((p, 0 , 10E9), (n, 0, 10E9), (m, 0, 10E9)) #when using infinity it goes weird, also integrating leads to higher thresholds than summation


a_mnp_denom = (((sp.sin(((m + 1) * sp.pi) / 2 * h)) ** 2) * ((sp.sin(((n + 1) * sp.pi) / W)) ** 2) * (
            (sp.sin(((p + 1) * sp.pi) / L)) ** 2) + K_hat)
a_mnp = A_mnp / a_mnp_denom  # Eqn. 18
U2 = sp.Function("U2")


U2 = a_mnp * ((sp.sin(((m + 1) * sp.pi) / 2 * h)) * x) * ((sp.sin(((n + 1) * sp.pi) / W)) * y) * (
            (sp.sin(((p + 1) * sp.pi) / L)) * z)  # Eqn. 13

x = np.arange(-h, h, 0.001)
y = np.arange(-h, h, 0.001)
z = np.arange(-h, h, 0.001)
f= sp.subs((U2), (x ,y ,z))

I currently get the error message: ValueError: subs accepts either 1 or 2 arguments. So that means I can't use the subs() method and replace() also doesn't work too well. Are there any other methods one can use?

Any help will be grateful, thank you!

1 Answers

Oscar is right: you are trying to deal with too much of the problem at once. That aside, Numpy and SymPy do not work like you think they do. What were you hoping to see when you replaced 3 variables, each with a range?

You cannot replace a SymPy variable/Symbol with a Numpy arange object, but you can replace a Symbol with a single value:

>>> from sympy.abc import x, y
>>> a = 1.0
>>> u = x + y + a
>>> u.subs(x, 1)
y + 2.0
>>> u.subs([(x,1), (y,2)])
4.0

You might iterate over the arange values, creating values of f and then doing something with each value:

f = lambda v: u.subs(dict(zip((x,y),v))) 
for xi in range(1,3): # replace range with your arange call
  for yi in range(-4,-2):
    fi = f((xi,yi))
    print(xi,yi,fi)

Be careful about iterating and using x or y as your loop variable, however, since that will then lose the assignment of the Symbol to that variable,

for x in range(2):
  print(u.subs(x, x)) # no change and x is no longer a Symbol, it is now an int
Related