How write the 'for' loop in mathematica

  • Context: Mathematica 
  • Thread starter Thread starter heinricht
  • Start date Start date
  • Tags Tags
    Loop Mathematica
Click For Summary

Discussion Overview

The discussion revolves around how to implement a 'for' loop in Mathematica, specifically translating a code snippet from another programming language into Mathematica syntax. Participants explore various methods to achieve this, including the use of functions, modules, and blocks. The conversation also touches on the use of the resulting data in differential equations.

Discussion Character

  • Technical explanation
  • Mathematical reasoning
  • Homework-related
  • Debate/contested

Main Points Raised

  • One participant presents a code snippet using a 'for' loop and asks for help translating it into Mathematica.
  • Another participant suggests using a 'Do' loop and provides an alternative implementation that creates alternating sequences of +g and -g.
  • A later reply proposes a more concise version of the code using 'Table' and 'Append' functions.
  • Participants discuss how to use the generated list 'p' as a function in an ODE, with one participant questioning the appropriateness of accessing list elements with non-integer indices.
  • Concerns are raised about the need for initial conditions in the ODE and the implications of using integer indices for accessing elements in the list 'p'.
  • Another participant suggests modifying the function to handle non-integer values of 't' by rounding down the index.
  • Further discussion includes the potential to redefine 'p' as a function rather than a list to avoid issues with indexing.
  • One participant shares a working program that uses 'p' as a constant and inquires about changing 'p' to a function in the context of an NDSolve problem.
  • There is a question about the time step in the NDSolve function and whether it can be adjusted.

Areas of Agreement / Disagreement

Participants express differing views on the best approach to implement the loop and how to handle the use of 'p' in subsequent calculations. There is no consensus on the optimal solution, and multiple competing methods are presented.

Contextual Notes

Participants highlight the importance of ensuring that the list 'p' is created before it is called in the module and the need to address the issue of non-integer indices when accessing elements from 'p'.

Who May Find This Useful

This discussion may be useful for individuals learning to program in Mathematica, particularly those transitioning from other programming languages and seeking to understand how to implement loops and manage data structures in Mathematica.

heinricht
Messages
7
Reaction score
0
I have a problem how write the 'for' loop in mathematica:

my previous code looks like:
count=0;
flag=0;
step=100;

For i=1:1000
count=count+1;
if count>step
if flag== 0;
flag=1;
else
flag=0;
end
count=1;
end
if flag==0
p(i)=g
else
p(i)=-g
end
end

I need to rewrite this in mathematica but I do not know how.
Is it possible to write as a function, module or block?
thanks
 
Physics news on Phys.org


Hi heinricht,

heinricht said:
I have a problem how write the 'for' loop in mathematica:

my previous code looks like:
count=0;
flag=0;
step=100;

For i=1:1000
count=count+1;
if count>step
if flag== 0;
flag=1;
else
flag=0;
end
count=1;
end
if flag==0
p(i)=g
else
p(i)=-g
end
end

I need to rewrite this in mathematica but I do not know how.
Is it possible to write as a function, module or block?
thanks

I think something like this would closely follow your program:

Code: 
count=0
flag=0
step=100

p=Table[0,{i,1,1000}]  (* add a semicolon after this to not display the zero matrix*)


Do[{
count=count+1,
If[ count > step , {If[ flag==0, flag=1, flag=0], count=1} ],
If[flag==0, p[[i]]=g, p[[i]]=-g]} , {i,1,1000}]

p   (* this is just displaying the result*)


It looks like your program is creating alternating sequences of +g and -g of 100 terms each. If that is your goal, you can more directly create it by using something like:

Code: 
pos=Table[g,{i,1,100}]

neg=Table[-g,{i,1,100}]

p={}

Do[ {p = Append[p,pos],
		p=Append[p,neg]},{i,1,5}]

p=Flatten[p]   (*removing nested brackets *)

p   (*display result*)

and you can easily shorten this quite a bit (but this form is understandable, I think).



As an example of a shorter program, you might try,

Code: 
p=g  Flatten[ 
      Table[ Append[ Table[1,{i,1,100}],Table[-1,{i,1,100}]],{i,1,5}]]

which does the same as the previous two programs.
 


Thanks a lot of. I have one more question:
I have equation: NDSolve[m (q[t] - M[t] - g (2 M[t] + 1) k) == M'[t], {t,0,1000}]
and I need to use p from previous program as a function of q[t]
I used p[[t]] which take t element from the list p and wrote
q[t_]=p[[t]]
is this correct?

or can I write sth like:

p[t_]:= Module[{step, avg, amp, freq, per, t1, t0, time},
amp = 0.1;
avg = 0.25;
freq = 10^8;
per = 1/freq;
t1 = 100/10^9;
t0 = 0;
time = (t1 - t0)/t;
step = per/time;

count = 0;
state = 0;
Do[{count = count + 1,
If[count > step, {If[state == 0, state = 1, state = 0],
count = 1}],
If[state == 0, p[] = avg + amp/2, p[] = avg - amp/2]}, {i,
1, t}];
p[[t]]];

and then use this function in ODE equation:
NDSolve[m (p[t] - M[t] - g (2 M[t] + 1) k) == M'[t], {t,0,1000}]
 
Last edited:


heinricht said:
Thanks a lot of. I have one more question:
I have equation: NDSolve[m (q[t] - M[t] - g (2 M[t] + 1) k) == M'[t], {t,0,1000}]

Don't you need an initial condition for this?

and I need to use p from previous program as a function of q[t]
I used p[[t]] which take t element from the list p and wrote
q[t_]=p[[t]]

I don't think this will work. As NDSolve is solving, you need to be able to get values for all t between 0 and 1000. However if you try to access q[1.5], for example, you'll get an error, since p[[t]] requires t to be an integer.

is this correct?

or can I write sth like:

p[t_]:= Module[{step, avg, amp, freq, per, t1, t0, time},
amp = 0.1;
avg = 0.25;
freq = 10^8;
per = 1/freq;
t1 = 100/10^9;
t0 = 0;
time = (t1 - t0)/t;
step = per/time;

count = 0;
state = 0;
Do[{count = count + 1,
If[count > step, {If[state == 0, state = 1, state = 0],
count = 1}],
If[state == 0, p[] = avg + amp/2, p[] = avg - amp/2]}, {i,
1, t}];
p[[t]]];

and then use this function in ODE equation:
NDSolve[m (p[t] - M[t] - g (2 M[t] + 1) k) == M'[t], {t,0,1000}]



Well, I'm not really certain what you are doing here, so I'll just make a few comments. You need to call the module something beside p, since p is the list inside the module (right now you have p[t_] on the left hand side of the equation, but also use p as the list called with p[]).

Also, you need to make sure the the list p is created before this module is called.

And here again there is the issue with non-integer values of t, which will cause an error with the module the way it is written. For example, let's say that you know that

p[[25]]=.2
p[[26]]=.3

Now when NDSolve gets to t=25.678, what value do you want it to call? There is no p[[256.678]]. Do you want it to be equal to p[[25]] or p[[26]], or something else?

For example, if you always want the index rounded down, you can change the last line of your module to this:


p[[ Floor[t] ]]];

which will just truncate the part of t after the decimal. Then when t=25.678, you would always get the value for p[[25]]. (Instead of Floor, you can use Round or Ceiling to get different behaviors.)



(Again, I'm not sure of your goal here, but you might want to consider not creating a list p, and instead just creating a function.)
 


my previuos program is below: I used p as a constant-number(0.25)

constant number: p,g0,ep,gamma,tau,m
number = 50000; (*number_steps*)

g10 = g0/(1 + ep F0);
int= NSolve[{(g10/g0 (2 M0 + 1) - 1) F0 == 0,
m (p - M0 - (g0 (2 M0 + 1)) F0) == 0} , {M0, F0}];
x0 = M0 /. int;
y0 = F0 /. int;
(*Initial conditions*)
m00 = x0[[3]]
f00 = y0[[3]]

g = g0/(1 + ep F[t]);
s = NDSolve [{(g/g0 (2 M[t] + 1) - 1) F[t] == F'[t],
m (p - M[t] - g (2 M[t] + 1) F[t]) == M'[t],
F[0] == f00, M[0] == m00},
{F[t], M[t]}, {t, 1, number}] ;
time = Table[{t/\[Gamma]}, {t, number}];
solF = Table[F[t] /. s, {t, number}];
solM = Table[M[t] /. s, {t, number}];
fo = Flatten[solF];
mo = Flatten[solM];
tg = Flatten[time];
final = Thread[{tg, fo}];
final2 = Thread[{tg, mo}];
ListLinePlot[final, PlotRange -> All]
ListLinePlot[final2, PlotRange -> All]

and this program is working fine and now I need use p in NDSolve as a function of p[t] using formula
amp = 0.1;
avg = 0.25;
f = 10^8;
per = 1/f;
t1 = 100/10^9;
t0 = 0;
time = (t1 - t0)/number;
step = per/time;

count = 0;
state = 0;
Do[{count = count + 1,
If[count > step, {If[state == 0, state = 1, state = 0],
count = 1}],
If[state == 0, p[] = avg + amp/2, p[] = avg - amp/2]}, {i,
1, t}]

beacuse p is also in NSolve where I have initial conditions from physical point of view I can use p0=0.01 so I do not have to calculate exactly value.

Can I change also in NDSolve time step=1? I gues that now the time step is 0.001 that I have funny thing p[[25.687]] if not could you help me with these function?
 


heinricht said:
my previuos program is below: I used p as a constant-number(0.25)

constant number: p,g0,ep,gamma,tau,m
number = 50000; (*number_steps*)

g10 = g0/(1 + ep F0);
int= NSolve[{(g10/g0 (2 M0 + 1) - 1) F0 == 0,
m (p - M0 - (g0 (2 M0 + 1)) F0) == 0} , {M0, F0}];
x0 = M0 /. int;
y0 = F0 /. int;
(*Initial conditions*)
m00 = x0[[3]]
f00 = y0[[3]]

g = g0/(1 + ep F[t]);
s = NDSolve [{(g/g0 (2 M[t] + 1) - 1) F[t] == F'[t],
m (p - M[t] - g (2 M[t] + 1) F[t]) == M'[t],
F[0] == f00, M[0] == m00},
{F[t], M[t]}, {t, 1, number}] ;
time = Table[{t/\[Gamma]}, {t, number}];
solF = Table[F[t] /. s, {t, number}];
solM = Table[M[t] /. s, {t, number}];
fo = Flatten[solF];
mo = Flatten[solM];
tg = Flatten[time];
final = Thread[{tg, fo}];
final2 = Thread[{tg, mo}];
ListLinePlot[final, PlotRange -> All]
ListLinePlot[final2, PlotRange -> All]

and this program is working fine and now I need use p in NDSolve as a function of p[t] using formula
amp = 0.1;
avg = 0.25;
f = 10^8;
per = 1/f;
t1 = 100/10^9;
t0 = 0;
time = (t1 - t0)/number;
step = per/time;

count = 0;
state = 0;
Do[{count = count + 1,
If[count > step, {If[state == 0, state = 1, state = 0],
count = 1}],
If[state == 0, p[] = avg + amp/2, p[] = avg - amp/2]}, {i,
1, t}]

beacuse p is also in NSolve where I have initial conditions from physical point of view I can use p0=0.01 so I do not have to calculate exactly value.

Can I change also in NDSolve time step=1?


I don't think you would want to do that; I believe mathematica needs to adjust the step size itself so as to get a good result.

I gues that now the time step is 0.001 that I have funny thing p[[25.687]] if not could you help me with these function?

The function needs to be a step function, right? There are many ways to do this, but I think a good way that never returns a zero would be something like (after you have already calculated step, avg, and amp):

Code: 
p[x_]:= If [ EvenQ [ Quotient [ x , step ] ] , avg + amp/2 , avg - amp/2 ]

Once that is defined, if I calculate:

Code: 
step=3;

avg=12;

amp=6;

Plot[p[x],{x,0,10} , PlotRange->{0,20}]

I get this:

http://img154.imageshack.us/img154/4494/squarewavedq2.jpg

So as long as you calculate step, avg, and amp before you run NDSolve I believe it will work.




(I also think you could use your module, but change the last line of the module to

p[[ Ceiling[t] ]] ];

to get the same behavior.)
 
Last edited by a moderator:


I used your function-see below- in my program. The step is fine (picture1) but when I used these function in NDSolve it does not work (picture 2): see the attachment


g0 = 1;
ep = 1/10^2;
Gamma = 0.98/10^-12;
Tau = 0.816/10^9;
m = 1/(Gamma*Tau);
number = 100 000;

g10 = g0/(1 + ep F0);
intensity = NSolve[{(g10/g0 (2 M0 + 1) - 1) F0 == 0,
m (avg - M0 - (g0 (2 M0 + 1)) F0) == 0} , {M0, F0}];
x0 = M0 /. intensity;
y0 = F0 /. intensity;
(*Initial conditions*)
Print["I0: ", m00 = x0[[3]]]
Print["N0: ", f00 = y0[[3]]]

amp = 0.1;
avg = 0.25;
freq = 10^8; (*Hz*)
per = 1/freq; (*period of modulated current*)
step = per/10^-12;

p[t_] := If[EvenQ[Quotient[t, step]], avg + amp/2, avg - amp/2];
Plot[p[t], {t, 0, 50000}, PlotRange -> {0, 0.4}]

g = g0/(1 + ep F[t]);
s = NDSolve [{(g/g0 (2 M[t] + 1) - 1) F[t] == F'[t],

m (p[t] - M[t] - g (2 M[t] + 1) F[t]) == M'[t],
F[0] == f00, M[0] == m00},
{F[t], M[t]}, {t, 1, number}] ;

Plot[F[t] /. s, {t, 0, number}, PlotRange -> All]
Plot[M[t] /. s, {t, 0, number}, PlotRange -> All]

I also used previous function: p[t_] := If[t > step, avg + amp/2, avg - amp/2];
but the result is only one step: see the picture 3 and 4
 

Attachments

  • data.JPG
    data.JPG
    15 KB · Views: 599
Last edited:


heinricht said:
I used your function-see below- in my program. The step is fine (picture1) but when I used these function in NDSolve it does not work (picture 2): see the attachment

What kind of result were you expecting? If you plot from 0 to 1000 (instead of from 0 to number) I believe you can see there are smooth oscillations; should there be something different?


(I noticed a few minor things, some of which I had to change to get your program to work for me. There was an extra space in the assignment for number (it was 100 000 instead of 100000), and I had to make the G in Gamma lowercase. Were those just typos?

Also, you had NDSolve go from 1 to number, and plotted the results from 0 to number. Was that on purpose?)
 


I have done these what I need exactly. I used simple function
p[t_] := If[t > step, If[t > 2*step, avg + amp/2, avg - amp/2], avg + amp/2];
and the result are below(data1,jpg)
Now I need to multiply steps only.

I tried used p[t_] as module and obtain the same result but it does not work. Can you tell me where might be the problem? Thanks.

p[t_] := Module[{step, avg, amp, freq, per, t1, t0, time, number,q},
amp = 0.1;
avg = 0.25;
freq = 10^8;
per = 1/freq;
t1 = 100/10^9;
t0 = 0;
number = 100000;
time = (t1 - t0)/number;
step = per/time;
count = 0;
state = 0;
Do[{count = count + 1,
If[count > step, {If[state == 0, state = 1, state = 0],
count = 1}],
If[state == 0, q[[t]] = avg + amp/2, q[[t]] = avg - amp/2]}, {i,
1, t}];
q[[Ceiling[t]]]
];
 

Attachments

  • data1.JPG
    data1.JPG
    8.5 KB · Views: 598
Last edited:
  • #10


heinricht said:
I have done these what I need exactly. I used simple function
p[t_] := If[t > step, If[t > 2*step, avg + amp/2, avg - amp/2], avg + amp/2];
and the result are below(data1,jpg)
Now I need to multiply steps only.

I tried used p[t_] as module and obtain the same result but it does not work. Can you tell me where might be the problem? Thanks.

p[t_] := Module[{step, avg, amp, freq, per, t1, t0, time, number,q},
amp = 0.1;
avg = 0.25;
freq = 10^8;
per = 1/freq;
t1 = 100/10^9;
t0 = 0;
number = 100000;
time = (t1 - t0)/number;
step = per/time;
count = 0;
state = 0;
Do[{count = count + 1,
If[count > step, {If[state == 0, state = 1, state = 0],
count = 1}],
If[state == 0, q[[t]] = avg + amp/2, q[[t]] = avg - amp/2]}, {i,
1, t}];
q[[Ceiling[t]]]
];


What is picture 2 from your last post? Is it a plot of the module p[t]? Or is it a plot of the results after you run NDSolve?


If it is just a plot of the module p[t], how did you plot it? Did you use

Plot[p[t],{t,0,50000}]

or something different? When I plot that for the module p, I get an identical result to your picture 1 (although there were a few errors about the list q that were reported, so I had to modify). It does not look to me like you need to build a list q inside the module, since only its last result is used.
 
  • #11


Hi in new year :),

the second picture is a result of NDSolve.
I did not plot the module p[t].
 
  • #12


Now I have a problem with number of steps in NDSolve.
My program:
avg = 0.5;
g10 = 1;
dg = 10^-3;
g20 = g10 - dg;
alpha = 2.6;
i = Sqrt[-1];
ep = 1/10^2;
Gamma = 0.98/10^-12;
Tau = 0.816/10^9;
m = 1/(Gamma*Tau);
beta = 2/3;
number = 30000;

amp = 0.2;
freq = 10^8;
tran = 10^-9;
per = 1/freq;
step = per/10^-12;
steptran = tran/10^-12;
phigh = avg + amp/2;
plow = avg - amp/2;
slope = (phigh - plow)/steptran;

c[t_] := If[t > step,
If[t > 2*step,
If[t < 3*step - steptran, phigh, phigh - (t - (3*step - steptran))slope ],
If[t < 2*step - steptran, plow, plow + (t - (2*step - steptran)) slope ] ],
If[t < step - steptran, phigh, phigh - (t - (1*step - steptran)) slope ] ];

Plot[c[t], {t, 0, number}, PlotRange -> {0.3, 0.65}]

Clear[M, F1, F2, t];
Mu = 1.0163;
Sigma = 0.75;
g17 = g10/(1 + ep F1[t] + ep beta F2[t]);
g27 = g20/(1 + ep beta F1[t] + ep F2[t]);
g1 = g17 1/(
1 - Sigma (1 - Mu)) (1 - Sigma] (2 M[t] + 1 - Mu));
g2 = g27 1/(
1 - Sigma (1 - Mu)) (1 + Sigma (2 M[t] + 1 - Mu));

s = NDSolve [{(g1/g10 (2 M[t] + 1) - 1) F1[t] == F1'[t],
(g2/g10 (2 M[t] + 1) - 1) F2[t] == F2'[t],
m (c[t] - M[t] - (2 M[t] + 1) (g1 F1[t]) + g2 F2[t]) == M'[t],
F2[0] == 0.01, F1[0] == 0.02, M[0] == 0.03},
{F1[t], F2[t], M[t]}, {t, 0, number}] ;
Plot[F1[t] /. s, {t, 0, number}, PlotRange -> All]
Plot[F2[t] /. s, {t, 0, number}, PlotRange -> All]

When I run I got message:
NDSolve::mxst: Maximum number of 10000 steps reached at the point t == \
1942.86059216402`. >>

is this connected with MaxSteps in NDSolve? or sth different
 
Last edited:

Similar threads

Mathematica Bar legend for a different function in mathematica
  • · Replies 1 ·
Replies
1
Views
3K
Mathematica Problems solving this differential equation for a Pendulum with Mathematica
  • · Replies 4 ·
Replies
4
Views
4K
MATLAB Conversion from Mathematica to Matlab
  • · Replies 1 ·
Replies
1
Views
2K
Mathematica How to check a particular solution of System of Linear ODEs?
  • · Replies 1 ·
Replies
1
Views
3K
Mathematica Plot gradient vector in ContourPlot
  • · Replies 3 ·
Replies
3
Views
4K
Mathematica Solve complex equation analytically
  • · Replies 19 ·
Replies
19
Views
3K
Mathematica Implementation of tessarines in Mathematica
  • · Replies 1 ·
Replies
1
Views
2K
Mathematica How to simplify expressions in Mathematica with log and exp in them?
  • · Replies 6 ·
Replies
6
Views
4K
Mathematica Numerical integration over a Green's function
  • · Replies 13 ·
Replies
13
Views
3K
Mathematica [Mathematica] How to use for loop to process a list of data
  • · Replies 6 ·
Replies
6
Views
5K