Mathematica problem: two nigh-identical codes

[MartinV]

Hello, everyone.

I have a rather bizzarre problem. I wrote a code that takes a matrix U[1] and calculates a matrix U[2] from it using simple math. I wrote three For loops, one to go across the U[n], the other two to read the matrix U[n] across all the lines and columns. I do it twice because it's a predictor/corrector kind of computaion. Pretty straightforward. Here's the code:

For[n = 1, n <= 10, n++, UU2 = ConstantArray[0, {dim, dim}];
UU3 = ConstantArray[0, {dim, dim}];
UU = U[n];
UU2[[1, All]] = UU[[1, All]]; UU2[[10, All]] = UU[[10, All]];
UU2[[All, 1]] = UU[[All, 1]]; UU2[[All, 10]] = UU[[All, 10]];
UU3 = UU2;
For[i = 2, i <= dim - 1, i++,
For[j = 2, j <= dim - 1, j++,
UU2[[i, j]] =
UU[[i, j]] - dt/dx*(G[UU[[i, j]]] - G[UU[[i - 1, j]]]) -
dt/dx*(H[UU[[i, j]]] - H[UU[[i, j - 1]]]) + dt*S[UU[[i, j]]];
]
]
For[i = 2, i <= dim - 1, i++,
For[j = 2, j <= dim - 1, j++,
UU3[[i, j]] =
0.5*(UU[[i, j]] + UU2[[i, j]] -
dt/dx*(G[UU2[[i + 1, j]]] - G[UU2[[i, j]]]) -
dt/dx*(H[UU2[[i, j + 1]]] - H[UU2[[i, j]]]) +
dt*S[UU2[[i, j]]]);
]
]
U[n + 1] = UU3;
]


The code above works. It gives me exactly what I need but it's suppose to be a bit unstable so I'm trying to write a different code now. The following code is nigh identical but for some reason Mathematica gives me all kinds of errors and fails to compute. Here's the faulty code:

For[n = 1, n <= 3, n++,
UU2 = ConstantArray[0, {dim, dim}];
UU3 = ConstantArray[0, {dim, dim}];
UU = U[n];
UU2[[1, All]] = UU[[1, All]]; UU2[[dim, All]] = UU[[dim, All]];
UU2[[All, 1]] = UU[[All, 1]]; UU2[[All, dim]] = UU[[All, dim]];
UU3 = UU2;
For[i = 2, i <= dim - 1, i++,
For[j = 2, j <= dim - 1, j++,
UU2[[i, j]] = UU[[i, j]] - dt/dx/2*(G[UU[[i, j]]] - G[UU[[i - 1, j]]]) + dt/2*S[UU[[i, j]]];
]
]
For[i = 2, i <= dim - 1, i++,
For[j = 2, j <= dim - 1, j++,
UU3[[i, j]] = 0.5*(UU[[i, j]] - UU2[[i, j]] - dt/2/dx*(G[UU2[[i + 1, j]]] - G[UU2[[i, j]]]) +
dt/2*S[UU2[[i, j]]]);
]
]
U[n + 1] = UU3;
]

The only difference I can see is that the calculation done in each step of the For loop is a bit shorter since it doesn't employ the H[UU2[[i,j]]] parts. I don't see how this could affect the overall scheme.

These are the errors Mathematica reports when I execute the second piece of code:

Set::write: Tag Times in Null Null U[2] is Protected. >>
Part::partw: Part All of U[2] does not exist. >>
Part::partd: Part specification U[2][[All,1]] is longer than depth of object. >>
Part::partd: Part specification U[2][[All,10]] is longer than depth of object. >>
Part::partw: Part 2 of U[2] does not exist. >>
Part::partw: Part 2 of U[2] does not exist. >>
General::stop: Further output of Part::partw will be suppressed during this calculation. >>
Part::partd: Part specification U[2][[1,2]] is longer than depth of object. >>
General::stop: Further output of Part::partd will be suppressed during this calculation. >>
Set::write: Tag Times in Null Null U[3] is Protected. >>
Set::write: Tag Times in Null Null U[4] is Protected. >>
General::stop: Further output of Set::write will be suppressed during this calculation. >>

Of course I clear all variables before running the second code and reload the definitions for G, H and S which are being used in the computation.

 

[Bill Simpson]

This
Set::write: Tag Times in Null Null U[2] is Protected.
probably means at least in part you are missing a couple of semicolons. Yes, sometimes Mathematica will let you get away with that, and sometimes it will bite you. Sometimes without a semicolon Mathematica thinks two things adjacent to each other, like i j, means you want to multiply i and j. But it can also think that For[...] x means you want to multiply the Null returned from For by x.

Stick semicolons between these two lines

For[i=2,i<=dim,i++...]]
For[i = 2, i <= dim - 1, i++,

and between these two lines

]
U[n + 1] = UU3;

that is at the bottom of your code.

That isn't going to fix all your problems and in particular not fix your "Part 2 of U[2] does not exist" problem, but it should get rid of all your "Tag Times in" problem.

Now you haven't shown the rest of your code where U[] is defined. Perhaps you mean that U is a matrix and not a function. If my guessing is correct then U[[]] instead of U[] may make some of your problems go away. Maybe U[] really is a function, but then U[n+1]=UU3 looks suspicious. I can't tell what you really meant.

 

[MartinV]

SEMICOLONS!

Blasted semicolons. That's all it took, at the ending of every ]. Thanks for the help, Bill. Sometimes I am such a ***.

 

[MartinV]

New problem, same code. Best not make another thread, I think.

I used the codes mentioned before in a predictor-corrector scheme. This calls the same method four times, each time feeding it the new result.

This is what the single method looks like:

For[i = 2, i <= dim - 1, i++,
For[j = 2, j <= dim - 1, j++,
(*Izračun prediktorja Lx1*)
UU2[[i, j]] =
UU[[i, j]] -
Dxy[UU[[i, j]]]*dt/dx/2*(G[UU[[i, j]]] - G[UU[[i - 1, j]]]) +
Dxy[UU[[i, j]]]*dt/2*S[UU[[i, j]]];
];
];
For[i = 2, i <= dim - 1, i++,
For[j = 2, j <= dim - 1, j++,
(*Izračun korektorja Lx1*)
UU3[[i, j]] =
0.5*(UU[[i, j]] - UU2[[i, j]] -
Dxy[UU[[i, j]]]*
dt/2/dx*(G[UU2[[i + 1, j]]] - G[UU2[[i, j]]]) +
Dxy[UU[[i, j]]]*dt/2*S[UU2[[i, j]]]);
];
];

I want to make that into a method or module (not sure how Mathematica calls it) so I can do something like this:
UU3 = method[UU, UU2];

(I used to use Matlab way more often than Mathematica so I'm still thinking in M-file ways.)

I tried using the := way of identifying my methods, but it only seems to work for one input parameter. The corrector part uses both UU and UU2 to get UU3. Somehow writing it method[UU, UU2]:=... doesn't work.

predCorr[UU_, UU2_] :=
UU2 = ConstantArray[0, {dim, dim}];
UU3 = ConstantArray[0, {dim, dim}];
UU2[[1, All]] = UU[[1, All]]; UU2[[dim, All]] = UU[[dim, All]];
UU2[[All, 1]] = UU[[All, 1]]; UU2[[All, dim]] = UU[[All, dim]];
UU3 = UU2;
For[i = 2, i <= dim - 1, i++,
For[j = 2, j <= dim - 1, j++,
(*Izračun prediktorja Lx1*)
UU2[[i, j]] =
UU[[i, j]] -
Dxy[UU[[i, j]]]*dt/dx/2*(G[UU[[i, j]]] - G[UU[[i - 1, j]]]) +
Dxy[UU[[i, j]]]*dt/2*S[UU[[i, j]]];
];
];
For[i = 2, i <= dim - 1, i++,
For[j = 2, j <= dim - 1, j++,
(*Izračun korektorja Lx1*)
UU3[[i, j]] =
0.5*(UU[[i, j]] - UU2[[i, j]] -
Dxy[UU[[i, j]]]*
dt/2/dx*(G[UU2[[i + 1, j]]] - G[UU2[[i, j]]]) +
Dxy[UU[[i, j]]]*dt/2*S[UU2[[i, j]]]);
];
];

I get all sorts of errors, particularly these:

Part::partd: Part specification UU[[1,All]] is longer than depth of object. >>
Set::noval: Symbol UU2 in part assignment does not have an immediate value. >>
Part::partd: Part specification UU[[10,All]] is longer than depth of object. >>
Set::noval: Symbol UU2 in part assignment does not have an immediate value. >>
Symbol::argx: Symbol called with 0 arguments; 1 argument is expected. >>
Set::noval: Symbol UU2 in part assignment does not have an immediate value. >>
Symbol::argx: Symbol called with 0 arguments; 1 argument is expected. >>
Set::noval: Symbol UU2 in part assignment does not have an immediate value. >>
Part::partd: Part specification UU[[2,2]] is longer than depth of object. >>
Part::partd: Part specification UU[[2,2]] is longer than depth of object. >>

I can tell that I can't simply copy/paste my working code into a delayed value scheme but I don't have a clue how to make it work.

I made sure all semicolons are there this time. ;)

 

[Bill Simpson]

Because of precedence I believe you are getting

predCorr[UU_, UU2_] := UU2 = ConstantArray[0, {dim, dim}];
(*and that is the end of your function definition and then *)
UU3 = ConstantArray[0, {dim, dim}];
(*and that is unrelated to your attempt to create a module *)

It is difficult for me to guess what you want given your description.
Perhaps this is something like what you are trying to do.

mymodule[UU_] := Module[{UU2, UU3},
For[i=2, i<=dim - 1, i++,
For[j=2, j<=dim-1, j++,(*Izračun prediktorja Lx1*)
UU2[[i,j]] = UU[[i,j]] - Dxy[UU[[i,j]]]* dt/dx/2* (G[UU[[i,j]]] - G[UU[[i-1,j]]]) + Dxy[UU[[i,j]]]*dt/2*S[UU[[i,j]]];];];
For[i=2, i<=dim-1, i++,
For[j=2, j<=dim-1, j++,(*Izračun korektorja Lx1*)
UU3[[i,j]]=0.5*(UU[[i,j]] - UU2[[i,j]] - Dxy[UU[[i,j]]]*dt/2/dx*(G[UU2[[i+1,j]]] - G[UU2[[i,j]]]) + Dxy[UU[[i,j]]]*dt/2*S[UU2[[i,j]]]);];];
UU3];
mymodule[someUU]

 

[Bill Simpson]

I made a mistake and left out the initialization of UU2 and UU3

mymodule[UU_] := Module[{UU2, UU3},
UU2 = ConstantArray[0, {dim, dim}]; (*Insert this *)
UU3 = ConstantArray[0, {dim, dim}]; (* and this *)
For[i=2, i<=dim - 1, i++,
...

Sorry

 

[MartinV]

No problem. I was just going to ask about that. ;)