Python code moving to Julia

New to Julia
1

Hey! So I have a Phyton code, that I want to run in Julia. I did all the formatting with function commands and things like that, but at a point, when I use an if cycle in a for cycle, it gives me the following error:
ERROR: LoadError: ArgumentError: indexed assignment with a single value to possibly many locations is not supported; perhaps use broadcasting .= instead?
I just couldn’t find a solution to this error. Could you help me with it, what is the problem and how can I change the code to make it work?

using LinearAlgebra
using PyPlot
using PyCall
np = pyimport("numpy")

# Parameters of inverted pendulum
M = 0.5     # mass of cart
m = 0.2     # mas of pendulum
L = 0.3     # length of pendulum
I = 0.006   # Moment of pendulum
b = 0.1     # friction coefficient
g = 9.80665 # gravity

# Approximated parameters
Ma = 0.15    # mass of cart
ma = 0.12     # mas of pendulum
La = 0.7     # length of pendulum
Ia = 0.016   # Moment of pendulum
ba = 0.15    # friction coefficient
g = 9.80665  # gravity

# State-space dimensions
nx = 4
nu = 1

# Define the nonlinear state space equations of inverted pendulum
function invertedPendulum(state, input_force)
    
    x, xdot, theta, theta_dot, F = state[0], state[1], state[2], state[3], input_force
    dxdt = np.zeros_like(state)
    D = (m+M)*(m*L^2+I)-m^2*L^2*np.cos(theta)^2
    dxdt[0] = xdot
    dxdt[1] = ((I+m*L^2)*(m*L*np.sin(theta)*theta_dot^2-b*xdot+F)+m^2*L^2*g*np.cos(theta)*np.sin(theta))/D
    dxdt[2] = theta_dot
    dxdt[3] = (-m*L*np.cos(theta)*(m*L*np.sin(theta)*theta_dot^2-b*xdot+F)-(m+M)*m*L*g*np.sin(theta))/D
    return dxdt
end

    # Approximated model to calculate the h(u)
function approx_model(t, td, tdd)
    
    #print('t, td, tdd, xd, xdd = ', t, td, tdd, xd, xdd)

    return hu
end

function RFPT_rule(u_approx, u, hu, K, B, A)
    #  need control parameters  K B A
    # If the norm of the error is greater then the limnit compute the deformation
    # (it is not near the Fixed Point)
    Bf = B*np.tanh(A*(hu-u_approx))
    out = (u+K)*(1+Bf)-K
    return out
end

error_limit = 1e-2
# Control Parameters For RFPT

K = 1e3
B = -0.5
A = 2e-3

Nt = 20000
dt = 0.001

TT = 2

x = np.zeros((nx, Nt))
u = np.zeros((nu, Nt))
time = np.arange(Nt)*dt

A1 = 0.7
w1 = np.pi / 8
theta_des = A1*np.sin(w1*time)
theta_dot_des = A1*w1*np.cos(w1*time)
theta_ddot_des = -A1*w1^2*np.sin(w1*time)

dxdt = np.zeros((nx,))
x0 = np.array([1, 0, 0.2, 0.3])
x[:,1] = x0

int_theta_error = 0
int_theta_dot_error = 0
upast = 0
hu = 0

for k in 1:(Nt-1)
    
    # Tracking Error Relaxation PI-type
    theta_error = theta_des[k]-x[2,k]
    theta_dot_error = theta_dot_des[k]-x[3,k]
    t = theta_des[k]
    td =  theta_dot_des[k]  + 2 * TT * theta_error     + TT^2 * int_theta_error 
    tdd = theta_ddot_des[k] + 2 * TT * theta_dot_error + TT^2 * theta_error

    # Calculate the input using approx model
    xdd = -((Ia+ma*La^2)*tdd+ma*g*La*np.sin(t))/(ma*La*np.cos(t))
    u_approx = (Ma+ma)*xdd + ma*La*tdd*np.cos(t)-ma*La*td^2*np.sin(t)

        if np.abs(upast-u_approx) < (1e-3)
             u[:,k] = u_approx
        else
              u[:,k] = RFPT_rule(u_approx, upast, hu, K, B, A)
        end
end

The error is in the last few line, which starts with the for cycle.

2

Assuming that u[:, k] is a slice of an array with multiple values and u_approx is a scalar, you need u[:, k] .= u_approx to set every entry of that slice to the scalar. In numpy, I think this automatically broadcasts (expands the scalar to the slice size).

2 Likes
3

Jules is likely right:

julia> x = rand(1, 1)
1×1 Matrix{Float64}:
 0.26754366641349914

julia> x[:, 1] = 0.5
ERROR: ArgumentError: indexed assignment with a single value to possibly many locations is not supported; perhaps use broadcasting `.=` instead?

(you don’t even need multiple values here, just slicing : on the left hand side)

Unrelatedly, why are you using all these numpy functions? As far as I can tell all this stuff is available in base Julia:

julia> zeros(3)
3-element Vector{Float64}:
 0.0
 0.0
 0.0

julia> sin(π)
0.0

julia> cos(π)
-1.0

julia> tanh(π)
0.99627207622075

julia> range(0, 1, 5)
0.0:0.25:1.0
2 Likes
4

Thanks, this seems to solve the issue.