Is Symbolics.jl incompatible with plotting, cannot convert Num to subtype, MethodError: no method matching Float64(::Num)

General Usage
, , ,
1

I was trying to write the analyical solution of a cantilever beam problem as
a symbolic expressions and I ran into a problem when I had just finished substituting the variables by numerical values and I wanted to plot these values.

I have substituted all my expressions completely with values, i.e, there are no more variables, however the problem seems to be due to symbolic expressions being of type Num and Plots/Makie/other packages, seem not to be able to convert Nums to floats, int, rationals and etc… And honestly neither can I,

I tried converting my values using float(), convert(), eval(), build_function(), I ofcourse looked up if Symbolics.jl had something build in, but I only came past the mention of symbolic_to_float() in “Symbolics.jl/test/utils.jl” file on GitHub, but I could not use/find this function as it does not seem to be listed anywere.

I could rewrite my code and avoid using symbolics, but I would like to know what the propper/recommended way of converting these Nums would be.

Below I have added some code to demonstrate the errors I get with the simplest code possible together with a test function to show that it is not a syntax error. How to reproduce

using Symbolics
@variables a b 
f(x) = a*x^2 + b
f_subs(x) = substitute(f(x), Dict(a=>1, b=>2))

f_subs(1), typeof(f_subs(1))
> (3, Num)

f_test(x) = x^2 + 2
f_test(1), typeof(f_test(1))
> (3, Int64)

Output with Plots.jl

using Plots
plot(f_test)
> # Plots normally

plot(f_subs)
> ERROR: MethodError: no method matching Float64(::Num)
>
> Closest candidates are:
  (::Type{T})(::Real, ::RoundingMode) where T<:AbstractFloat
   @ Base rounding.jl:207
  (::Type{T})(::T) where T<:Number
   @ Core boot.jl:792
  (::Type{T})(::SymbolicUtils.Symbolic) where T<:Union{AbstractFloat, Integer, Complex{<:AbstractFloat}, Complex{<:Integer}}
   @ Symbolics C:\Users\userName\.julia\packages\Symbolics\VIBnK\src\Symbolics.jl:146
  ...

Output with GLMakie.jl (Plots is not loaded)

using GLMakie
lines(f_test.(-5:5)) 
> # Plots normally

lines(f_subs.(-5:5))
> `Makie.convert_arguments` for the plot type Lines and its conversion trait PointBased() was unsuccessful.

 The signature that could not be converted was:
 ::Vector{Num}

 Makie needs to convert all plot input arguments to types that can be consumed by the backends (typically Arrays with Float32 elements).
 You can define a method for `Makie.convert_arguments` (a type recipe) for these types or their supertypes to make this set of arguments convertible (See http://docs.makie.org/stable/documentation/recipes/index.html).

 Alternatively, you can define `Makie.convert_single_argument` for single arguments which have types that are unknown to Makie but which can be converted to known types and fed back to the conversion pipeline.

 Substrace:
 ...

 caused by: MethodError: no method matching Float32(::Num)
 ...
 etc.

```
2

The problem here seems to be that the plotting variable is not a symbolic variable. The plot recipe for Num converts the Num to a function of its variables, but f_subs has none.

Is there a reason you’ve defined f as a numeric function with two symbolic variables rather than a purely symbolic function?

1 Like
3

I did that because I believe it is currently not possible to plot symbolic function in terms of variables, and besides that I would like to be able to substitute different beam dimensions and loading forces, to compare it to fem results (which I can if I substitute the dims and x-coordinate), but I also wanted to use the equations obtained to plot the deflection/angle/moment/shearforce diagram.

I technically speaking defined x (coordinate direction of the beam) also as a variable in order to read the short simplified version of the equation, but whether it is substituted() or passed as a function(x) varible, it does not seem to matter for plotting, since even if all the values I wanted to plot would be, e.g., Float64 (they are all Nums, I say this to demonstrate the issue as I perceived), if one of the values I pass as an column Vector is a Num(subbed), there seems to be no method to convert the Num into a subtype that is usable by the functions or use it as is

4

I did that because I believe it is currently not possible to plot symbolic function in terms of variables

@variables a b x
f = a*x^2+b
f_subs = substitute(f, Dict(a=>1, b=>2))
plot(f_subs)

works fine on my machine. Is this not workable for your use case?

5

Here’s another solution with minimal change of the original code. Just use

plot(Symbolics.value ∘ f_subs)

and a plot shows up. This is because the function Symbolics.value converts a number of the Num type (defined by Symbolics) to Julia’s native number type.

1 Like
6

My apologies, what I meant by it is not possible is that you cannot get the x/y-axis to (easily) include ticks like qL^4/(8EI) or something like that and not that it is not possible to use the variable x.

With that being said: Symbolics.value() seems to be exactly what I was looking for, so thank you very much. It took me a bit to get used to the syntax, but I managed to figure it out thankfully. I am not fully happy with how I had to use two different functions for each V/M/phi/u, since it will likely become a large mess if I tried to add different loading conditions. For now, my implementation is sufficient as I do not have infinite time, however I have added it below in case anyone is curious or would like to improve it.

using Symbolics

@variables L q d b #beam length, distributed load, height/thickness/y, width/x
@variables E #Young's modulus

x_load_0 = L*3//7
q_length = L-x_load_0 #L*4//7

A = d*b #height*width
Ixx = 1//12*d^3*b
#Iyy = 1//12*d*b^3


Ry_0 = q*q_length
Mz_0 = -q*q_length*(x_load_0 + 1//2*q_length)

# dV/dx = q
# dM/dx = V/(E*Ixx)
# dphi/dx = M/(E*Ixx)
# du/dx = phi/(E*Ixx)
@variables x

dVdx = [0,-q] # x = [0, x_load_0], [x_load_0, L]

V = [Ry_0, q*(L-x)]/(E*Ixx) # x = [0, x_load_0], [x_load_0, L]

C1 =  Mz_0 + Ry_0*x_load_0 - q*(L*x_load_0-x_load_0^2//2) #From solving +C with I.C. M(x_load_0)=Mz_0+Ry_0*x_load_0 or (when everything done right) equivalently M(L)=0
M = [(Mz_0 + Ry_0*x)
     (C1 + q*(L*x-x^2//2))]/(E*Ixx) # x = [0, x_load_0], [x_load_0, L]

C2 = Mz_0*x_load_0 + Ry_0*x_load_0^2//2 - q*(L*x_load_0^2//2 - x_load_0^3//6) - C1*x_load_0
phi = [(Mz_0*x + Ry_0*x^2//2)
       (C2 + C1*x + q*(L*x^2//2 - x^3//6))]/(E*Ixx)

C3 = Mz_0*x_load_0^2//2 + Ry_0*x_load_0^3//6 - q*(L*x_load_0^3//6 - x_load_0^4//24) - C1*x_load_0^2//2 - C2*x_load_0
u = [(Mz_0*x^2//2 + Ry_0*x^3//6)
     (C3 + C2*x + C1*x^2//2 + q*(L*x^3//6 - x^4//24))]/(E*Ixx)


n = 100; # number of elements per side
vars = [L,q,d,b,E]
vals = [10, 4, 1//10, 1//10, 200e9] #L,q,d,b,E
subs_Dict = Dict(vars .=> vals)
subs(v_) = substitute(v_, subs_Dict)
subs(v_,x_) = substitute(v_, subs_Dict)
subs(v_,x_::Vector) = [substitute(v_, Dict(subs_Dict..., Dict(x => i)...)) for i in x_]

# Plot
using Plots

V_subs1(x_) = Symbolics.value(substitute(subs(V[1]), Dict(x=>x_)))
V_subs2(x_) = Symbolics.value(substitute(subs(V[2]), Dict(x=>x_)))

M_subs1(x_) = Symbolics.value(substitute(subs(M[1]), Dict(x=>x_)))
M_subs2(x_) = Symbolics.value(substitute(subs(M[2]), Dict(x=>x_)))

phi_subs1(x_) = Symbolics.value(substitute(subs(phi[1]), Dict(x=>x_)))
phi_subs2(x_) = Symbolics.value(substitute(subs(phi[2]), Dict(x=>x_)))

u_subs1(x_) = Symbolics.value(substitute(subs(u[1]), Dict(x=>x_)))
u_subs2(x_) = Symbolics.value(substitute(subs(u[2]), Dict(x=>x_)))

x_load_0_rat = Symbolics.value(subs(x_load_0))
x_step_L = collect(0:x_load_0_rat*2/n:x_load_0_rat);
x_step_R = collect(x_load_0_rat:(subs_Dict[L]-x_load_0_rat)*2/n:subs_Dict[L]);

PV1 = plot(x_step_L, V_subs1.(x_step_L));
plot!(PV1, x_step_R, V_subs2.(x_step_R));

PM1 = plot(x_step_L, M_subs1.(x_step_L));
plot!(PM1, x_step_R, M_subs2.(x_step_R));

Pphi1 = plot(x_step_L, phi_subs1.(x_step_L));
plot!(Pphi1, x_step_R, phi_subs2.(x_step_R));

Pu1 = plot(x_step_L, u_subs1.(x_step_L));
plot!(Pu1, x_step_R, u_subs2.(x_step_R));

titles = ["V shear force" "M moment" "phi angle" "u displacement"]
fig = Plots.plot(PV1, PM1, Pphi1, Pu1, xlims=(0, subs_Dict[L]*1.1), layout=4, title=titles)