Converting a Matlab code into Julia

New to Julia
1

Hello there, i’m trying to convert a Matlab code into Julia, but no success so far. Here it is. I’m loading some numeric data from a txt file.

workspace()

dados = readdlm("C:/Users/Muril/AppData/Local/atom/app-1.28.0/Arquivos Feitos/Dados.txt")

Nt = dados[2,12]; 
maxiter = dados[2,13]; 
tol = dados[2,14]; 
U0 = dados[2,15]*1e3; 
delta = dados[2,16]; 
de = dados[2:Nt+1,1]; 
para = dados[2:Nt+1,2]; 
Lcomp = dados[2:Nt+1,3];
r = dados[2:Nt+1,4]; 
x = dados[2:Nt+1,5]; 
R = (r.*Lcomp); 
X = (x.*Lcomp); 
PLlin = (dados[2:Nt+1,7])*1e3;
QLlin = (dados[2:Nt+1,8])*1e3; 
PLnlin = (dados[2:Nt+1,9])*1e3;
QLnlin = (dados[2:Nt+1,10])*1e3; 
Qshunt = (dados[2:Nt+1,11])*1e3; 
Perc = 105; 
Qmod = dados[2,18]*1e3; 
kj = dados[2,19]; 
fq = dados[2,20]; 
hmax = 25; 
f = 60; 

PLtotal = (PLlin + PLnlin); 
QLtotal = (QLlin + QLnlin - Qshunt); 



maiornivel = 1; 
nivel = zeros(1,Nt); 
nivel[find(de .== 0)] = maiornivel;
i = 0;
while nivel[sum(nivel .== 0)] >= 0
for i = 1:Nt
if nivel[i] == maiornivel
nivel[find(de .== i)] = maiornivel+1;
end
end
maiornivel = maiornivel+1;
if nivel != 0
    break
end
end

Here, Nt = 17, and my “nivel” vector should be (got this from Matlab):

nivel = [1 2 3 4 5 6 7 8 9 4 3 4 5 5 6 6 7]

But when i run the above code in Juno, i get:

nivel = [1 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]

So, in my while loop, only one iteration is being made, but i don’t know why. Any help would be trully appreciated, thanks in advance!

2

Not sure what de is and therefore what exactly is in nivel when you start, but it seems to me that the if nivel != 0 might break the loop prematurely?

Have a look at the @show macro which might help you figure out what’s happening inside the loop

3

Comparison between a vector and a scalar works differently in Matlab and Julia. You probably want to replace

if nivel != 0

with

if all(nivel .!= 0)

Alternatively something similar if I didn’t get the Matlab semantics entirely right.

1 Like
4

Taking this a little deeper. A vector and a scalar are never equal in Julia, even if the vector has only one element and it has the same value as the scalar.

julia> [0] != 0
true

Therefore your loop is always terminated after one iteration. However, what really happens in Matlab is that the vector/scalar comparison is automatically vectorized. Julia can also vectorize the comparison, but you have to be explicit about it with a dotted operator.

julia> [0, 1, 2, 0] .!= 0
4-element BitArray{1}:
 false
  true
  true
 false

Now comes the second difference between Julia and Matlab. If you try to use this result in an if statement

julia> if [0, 1, 2, 0] .!= 0
       end
ERROR: TypeError: non-boolean (BitArray{1}) used in boolean context

you find that Julia complains, because it requires the expression to yield a (scalar) boolean, whereas Matlab happily accepts anything and decides whether it is to be considered as false or true. Now it’s not terribly obvious whether e.g. a vector of booleans should be considered true if some element is true or it should require all elements to be true. Matlab has a rule about that, which you need to remember if you use this kind of expression. Julia instead requires you to be explicit about it, using any or all to reduce the vector to a scalar.

julia> any([0, 1, 2, 0] .!= 0)
true

julia> all([0, 1, 2, 0] .!= 0)
false

julia> all([1, 1, 2, 3] .!= 0)
true

Now we can use this in an if statement, which leads to the if all(nivel .!= 0) suggestion.

12 Likes
5

Hello again! Thanks everyone for the fast responses! The command

all(nivel .! = 0)

did the work!

But i’m sorry, i’m so confused and i got stuck again in another part of my code… if anyone could help me out. Here it is:

for iter = 1:maxiter

for i = ordem
k = find(de .== i);

P[i] = PLtotal[i]+sum(P[k]+DP[k]);

end
end

Here the vector “ordem” is this 1x17:

ordem = [9 8 7 17 6 15 16 5 13 14 4 10 12 3 11 2 1]

and the vector “de” is this 1x17:

de = [0 1 2 3 4 5 6 7 8 3 2 11 12 12 14 14 16]

So, in the first iteration of the inner for, i = 9.

“PLtotal”, “P” and “DP” are 1x17 vectors with some numeric data.

In MATLAB, when i get into the line “k = find(de .== i)”, i get an empty matrix 0x0 (because theres no 9 in the “de” vector). Looking at the definition of the “sum” function in Matlab, it says this:

If A is an empty 0-by-0 matrix, then sum(A) returns 0

So, thats why it does that

sum(P[k]+DP[k])

just fine in there, returning me 0.

But in Julia, i get an index error, because k isn’t defined, since the expression “find(de .== i)” gives me a 0-element array.

Is there anyway i can fix this? I mean, convert the empty array into 0 just like Matlab does? Or any other viable solution.

7

Note that this allocates a temporary array. You can also do all(!iszero, nivel) or simply !iszero(nivel).

(In general, Matlab-style code that tries to cram everything into vector operations with no regard for allocation of temporary arrays or the number of loops is unlikely to take much if any advantage of Julia’s speed potential.
Not to mention the Matlab idiom of long scripts operating on globals. I’m guessing that once you get your code working you’ll find that it is no faster than Matlab, and maybe slower… at that point come back and people can help with performance tips.)

2 Likes
9

I’m very unclear on what the problem is. If I understand correctly you say that something along these lines gives an error, which I can’t replicate.

julia> de = [0 1 2 3 4 5 6 7 8 3 2 11 12 12 14 14 16]
1×17 Array{Int64,2}:
 0  1  2  3  4  5  6  7  8  3  2  11  12  12  14  14  16

julia> P = zeros(Int,1,17)
1×17 Array{Int64,2}:
 0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0

julia> DP = zeros(Int,1,17)
1×17 Array{Int64,2}:
 0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0

julia> i = 9
9

julia> k = find(de .== i)
0-element Array{Int64,1}

julia> sum(P[k]+DP[k])
0

It’s always helpful to give a complete but as small example as possible of what goes wrong so people can just paste it into the REPL and experiment with it themselves.

10

I got it, it was because my “ordem” vector had float numbers in it. I had to convert them to int, so they could be used in the index of the other vectors.

Ok, almost done… but now i got some problems with return of functions. Heres the function i converted from Matlab to Julia.

function MSP(Nt,maxiter,tol,U0,delta,de,para,Lcomp,r,x,PLlin,PLnlin,QLlin,QLnlin,Qshunt,Perc)

R = (r.*Lcomp)';
X = (x.*Lcomp)';
PLtotal = (PLlin + PLnlin);
QLtotal = (QLlin + QLnlin - Qshunt);

maiornivel = 1;
nivel = zeros(1,Nt);
nivel[find(de .== 0)] = maiornivel;
i = 0;
while nivel[sum(nivel .== 0)] >= 0
for i = 1:Nt
if nivel[i] == maiornivel
nivel[find(de .== i)] = maiornivel+1;
end
end
maiornivel = maiornivel+1;
if all(nivel .!= 0)
    break
end
end

ordematual = 1; 
ordem = zeros(1,Nt);
for i = maiornivel:-1:1
for m = 1:Nt
if nivel[m] == i
    ordem[ordematual] = m;
ordematual = ordematual+1;
end
end
end

DP = zeros(1,Nt);
DQ = zeros(1,Nt);
DPTotalAnt = 0;
DQTotalAnt = 0;
P = zeros(1,Nt);
Q = zeros(1,Nt);
u = zeros(1,Nt);
gama = zeros(1,Nt);
DPTotal = 0;
DQTotal = 0;

ordem = round.(Int,ordem); 

for iter = 1:maxiter

for i = ordem

k = find(de .== i);
P[i] = PLtotal[i]+sum(P[k]+DP[k]);
Q[i] = QLtotal[i]+sum(Q[k]+DQ[k]);

end

i = ordem[Nt];
A = R[i]*P[i]+X[i]*Q[i]-U0^2/2;
C = (R[i]^2+X[i]^2)*(P[i]^2+Q[i]^2);
B = sqrt(A^2-C);
u[i] = (Perc/100)*sqrt(B-A);
gama[i] = delta - asin((P[i]*X[i]-Q[i]*R[i])/(U0*u[i]));

for i = 2:Nt
A = R[i]*P[i]+X[i]*Q[i]-u[de[i]]^2/2;
C = (R[i]^2+X[i]^2)*(P[i]^2+Q[i]^2);
B = sqrt(A^2-C);
u[i] = sqrt(B-A);
gama[i] = gama[de[i]] - asin((P[i]*X[i]-Q[i]*R[i])/(u[i]*u[de[i]]));
end

DP = R.*(P.^2+Q.^2)./(u.^2);
DQ = X.*(P.^2+Q.^2)./(u.^2);
DPTotal = sum(DP);
DQTotal = sum(DQ);

variacao = abs(DPTotal-DPTotalAnt);
DPTotalAnt = DPTotal;
DQTotalAnt = DQTotal;
if variacao < tol
break
end
end

end

In Matlab, i could just do:

function [u,gama,P,Q,DP,DQ,DPTotal,DQTotal] = MSP(Nt,maxiter,tol,U0,delta,de,para,Lcomp,r,x,PLlin,PLnlin,QLlin,QLnlin,Qshunt,Perc)

And i would get the vectors u,gama,P,Q,DP,DQ,DPTotal,DQTotal. How can i get those vectores in Julia? I tried return u,gama,P,Q,DP,DQ,DPTotal,DQTotal just before the last end of the function, but it doesn’t seem to work that way.

11

That ought to work—there don’t appear to be any other exit points of this function.

12

Guessing wildly your problem might be using Matlab syntax to retrieve the multiple outputs. At any rate, this is how returning multiple values typically looks in Julia.

julia> function f()
           return 1, 2
       end
f (generic function with 1 method)

julia> a, b = f()
(1, 2)

julia> a
1

julia> b
2
13

Thansk a lot Gunnar! This worked for me!

14

The replies here are obviously more immediately useful but as you continue in converting MATLAB code to Julia, if you’d like to document what you learn, I’ve created a wiki geared towards this that you could contribute to. See MATLAB to Julia wiki.

1 Like
15

Hello nooblia, very nice iniciative!! As soon as i end my code, i’ll look into it, and share everything i’m learning here in the discourse forums.

So guys, i have another question now… here’s a part of my code:

for i = 1:17
    i4[i] = conj((PLnlin[i]+im*QLnlin[i])/u[i]);
end

iaux1 = i4;

C = [100 0 0 0 20 0 14.3 0 0 0 9.1 0 7.7 0 0 0 5.9 0 5.3 0 0 0 4.3 0 4];
Ih = [];
for h = 3:2:25
    for i = 1:17
        if PLnlin[i] != 0
          i4[i] = (C[h]/100)*i4[i];
        end
    end
    Ih = [Ih;i4];
    i4 = iaux1;
end

And here we got:

PLnlin = [0.0 0.0 0.0 0.0 0.0 3.0e6 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0]

So,

PLnlin[6] = 3e6

But when i get into i = 6, that if isn’t working. I’ve tried the following:

julia> PLnlin[6] != 0
true

Why it’s not working when i run my code? If you need any other info to analyze this, just let me know. Btw, QLnlin and u are also 1x17 vectors with numeric data.

.

16

Being able to copy paste the code and have it run helps.

iaux1 = i4;

looks odd, if you want a copy you need iaux1 = copy(i4).

17

Hello Kristoffer!

That copy(i4) worked! That’s pretty weird for me, that using iaux1 = i4; doesn’t work. Guess i’m really used to Matlab…

Thanks for your help!

18 
A = [1,2,3]

means you have a variable A bound to an object [1,2,3]:

A --------> [1,2,3]

When you write B = A you get

A ----------v
            [1,2,3]
B-----------^

So if we do e.g. A[1] = 2 we will have

A ----------v
            [2,2,3]
B-----------^

and thus B[1] == 2.

If you do B = copy(A) you get

A --------> [1,2,3]
B --------> [1,2,3]

Now if you do A[1] = 2 you get

A --------> [2,2,3]
B --------> [1,2,3]

so in this case B[1] = 1.

28 Likes
19

These ASCII diagrams should go in a FAQ or somewhere.

13 Likes
20

Hello again!

Supose i have the following vector and a maximum value of 5 (to compare with the values in the vector).

v = [1 2 10 12 3 4 5]

How can i print something to show in wich positions of the vector v, has values that are higher than 5? I was thinking something like:

if any(v .> 5)
   println("detect wich values are higher than 5 and print their position")
else
  println("no values exceded the maximum limit")
end

but i’m not sure whats suposed to go in the argument of the println function.

21 
julia> using Random # needed for Shuffle on Julia 0.7

julia> x = shuffle!(collect(1:10));

julia> x'
1×10 LinearAlgebra.Adjoint{Int64,Array{Int64,1}}:
 3  8  6  7  5  1  4  10  2  9

julia> [x[i] for i ∈ 1:10 if x[i] > 5]
5-element Array{Int64,1}:
  8
  6
  7
 10
  9

julia> [i for i ∈ 1:10 if x[i] > 5]
5-element Array{Int64,1}:
  2
  3
  4
  8
 10

So,

julia> y = [i for i ∈ 1:10 if x[i] > 5];

julia> if length(y) > 0
           println("These indices correspond to values greater than 5:\n", y)
       else
           println("no values exceded the maximum limit")
       end
These indices correspond to values greater than 5:
[2, 3, 4, 8, 10]

If you’d rather not allocate the vector y, you could just use a for loop and print. Note that any(v .> 5) allocates, anyway.

2 Likes
22

Great! Thanks a lot Elrod!