RE: Re: Re: Particular structure 2
- To: mathgroup at smc.vnet.net
- Subject: [mg33733] RE: [mg33713] Re: [mg33627] Re: Particular structure 2
- From: "Wolf, Hartmut" <Hartmut.Wolf at t-systems.com>
- Date: Wed, 10 Apr 2002 00:49:18 -0400 (EDT)
- Sender: owner-wri-mathgroup at wolfram.com
> -----Original Message----- > From: Yas [mailto:yast at optushome.com.au] To: mathgroup at smc.vnet.net > Sent: Tuesday, April 09, 2002 7:03 AM > Subject: [mg33733] [mg33713] [mg33713] Re: [mg33627] Re: Particular structure 2 > > > Dear Allan and Hartmut, > Thank you very much for your replies to my questions and > the solutions > you have provided. As has been pointed out, a best solution is hardly > feasible without further details on the problem. The last reply from > Allan ([mg33656]) containing the construct; > > x= y=Table[1,{400}];f=Plus;F= Times; > > Table[ReleaseHold[F@@#],{j, Length[y]}]&[ > Table[f[x[[i]],Hold[y[[j]]],z],{i,Length[x]}]];//Timing > > is again another factor of two improvement in terms of speed but > slightly less efficient than Hartmut Wolfs Case 6 (mg33657), > > Plus @@ Map[Function[{yy}, Evaluate[F[Map[Function[{x}, f[x, yy, z]], > x]]]], y] // Timing > > in terms of memory use. > > I was particularly encouraged by the speed of Thread, but > unfortunately > the structure of the resulting matrix gets jumbled up. I > couldn't see a > way around it. Transpose restructures it in a way where the > 3x3 and 6x6 > relations get messed up. > > In my problem (the simplest case) f is the equivalent of a > 3x3 matrix > and is defined as, evalmat[x_,y_,z_] = mat[x, y, z], and mat has been > built up from matrix elements containing functions depending > on x, y and > z. For fixed values of y and z and with x, a list of 500 or more > elements, I use Fold or FoldList to construct a final vector > or a list > of vectors. So the new vector, newvec = Fold[#2.#1&, vec, > g1], where vec > ={0,0,1} usually, and g1 is a list built by mapping x onto > evalmat (g1 = > evalmat[#, y, z]&/@x). Incidentally I found using Outer to be > slightly > quicker for mapping x onto evalmat with y held fixed which has had me > stumped as to why. > > For the same list of x values and with y now variable as well my > natural inclination was to create a large array of g1's with > incrementing y which can be done by adapting either of the two > constructs above, > > e.g. g1= Map[Function[{yy}, Evaluate[Map[Function[{x}, > evalmat[x, yy, > z]], x]]], y] > or g1 = Table[ReleaseHold[#],{j, Length[y]}]&[ > Table[evalmat[x[[i]],Hold[y[[j]]],z],{i,Length[x]}]]; > > The building and evaluation of g1 is the slowest part if y is > variable > as well, and for sufficient resolution where any spurious jumps are > avoided the number of increments of y has to be large, usually > 256 > points. But once g1 is built, and because the Fold part is fast the > construct below is not limiting. > > profile = Table[Fold[#2.#1&, vec, g1[[k]]], {k, > Length[y]}]; > > The elements of the vector can be extracted and plotted. > > Now that all this has been achieved, I immediately > applied it to a > larger matrix (6x6) only to find that I am back to square one > in terms > of time and memory. In the 3x3 case all the matrix elements have set > delayed assignment and the functions inside the elements of > the matrix > also have delayed assignment. The matrix mat is also delayed but I > figure using evalmat, rather than mat overcomes re-evaluation > for each > new x, y or z. Replacing := with = for the matrix element definitions > and mat so that they are evaluated immediately did not seem to make a > difference to the evaluation times, so what else can be done > to improve > the time for evaluation. > > > Regards > > Yas > Yas, after all, best would be to make a clear cut, communicable model of your problem. I'm not shure, but try: evalmat[x_, y_, z_] := {{x^2, x y, x z}, {0, y^2, y z}, {0, 0, z^2}} (just to be concrete here, it's neccessary, to get at the problem) xx = {x1, x2, x3}; yy = {y1, y2}; (symbolic, and short to verify correctness) Now if this is a solution to your problem (my best guess of your description): Function[y, Fold[#2.#1 &, {0, 0, 1}, evalmat[#, y, z] & /@ xx]] /@ yy {{x3*z^5 + x3^2*(x1*x2^2*z + x2*y1^2*z + x2*z^3) + x3*y1*(y1^3*z + y1*z^3), y1*z^5 + y1^2*(y1^3*z + y1*z^3), z^6}, {x3*z^5 + x3^2*(x1*x2^2*z + x2*y2^2*z + x2*z^3) + x3*y2*(y2^3*z + y2*z^3), y2*z^5 + y2^2*(y2^3*z + y2*z^3), z^6}} What is done here? First the Function is defined, but the body is not yet evaluated. Then each value y of list yy ist inserted into the body, then body is evaluated. Thus whenever evalmat is reached, all its arguments are numeric (when you do it within your application). This is probably most memory efficient. BTW instead of Fold[#2.#1 &, {0, 0, 1}, evalmat[#, y, z] & /@ xx] {x3*z^5 + x3^2*(x1*x2^2*z + x2*y^2*z + x2*z^3) + x3*y*(y^3*z + y*z^3), y*z^5 + y^2*(y^3*z + y*z^3), z^6} you may do in a single scan Fold[evalmat[#2, y, z].#1 &, {0, 0, 1}, xx] {x3*z^5 + x3^2*(x1*x2^2*z + x2*y^2*z + x2*z^3) + x3*y*(y^3*z + y*z^3), y*z^5 + y^2*(y^3*z + y*z^3), z^6} If you compare the double map with Function[y, Evaluate[Fold[#2.#1 &, {0, 0, 1}, evalmat[#, y, z] & /@ xx]]] /@ yy {{x3*z^5 + x3^2*(x1*x2^2*z + x2*y1^2*z + x2*z^3) + x3*y1*(y1^3*z + y1*z^3), y1*z^5 + y1^2*(y1^3*z + y1*z^3), z^6}, {x3*z^5 + x3^2*(x1*x2^2*z + x2*y2^2*z + x2*z^3) + x3*y2*(y2^3*z + y2*z^3), y2*z^5 + y2^2*(y2^3*z + y2*z^3), z^6}} Then the body of the outer function is evaluated first with *symbolical* y (which may burden the memory compared with pure numerical evaluation), but the inner Map and Fold is evaluated only once (for all values y of yy). This normally increases speed. This is quite similar to Fold[#2.#1 &, {0, 0, 1}, evalmat[#, y, z] & /@ xx] /. List /@ Thread[y -> yy] {{x3*z^5 + x3^2*(x1*x2^2*z + x2*y1^2*z + x2*z^3) + x3*y1*(y1^3*z + y1*z^3), y1*z^5 + y1^2*(y1^3*z + y1*z^3), z^6}, {x3*z^5 + x3^2*(x1*x2^2*z + x2*y2^2*z + x2*z^3) + x3*y2*(y2^3*z + y2*z^3), y2*z^5 + y2^2*(y2^3*z + y2*z^3), z^6}} Now why doesn't Thread work in your case? Regard: Evaluate[Fold[#2.#1 &, {0, 0, 1}, Thread[evalmat[xx, #, z]]]] & /@ yy Dot::"rect": "Non-rectangular tensor encountered." Dot::"rect": "Non-rectangular tensor encountered." Dot::"rect": "Non-rectangular tensor encountered." General::"stop": "Further output of \!\(Dot :: \"rect\"\) will be suppressed \ during this calculation." {{{x1*z, x2*z, x3*z}, y1*z, z^2} . {{x1*y1, x2*y1, x3*y1}, y1^2, 0} . {{x1^2, x2^2, x3^2}, 0, 0} . {0, 0, 1}, {{x1*z, x2*z, x3*z}, y2*z, z^2} . {{x1*y2, x2*y2, x3*y2}, y2^2, 0} . {{x1^2, x2^2, x3^2}, 0, 0} . {0, 0, 1}} The error comes because now evalmat is evaluated *before* Thread becomes effective. **This is a very common error when using Thread**, quite often discussed in MathGroup. What we have to do is to hold evalmat by wrapping it into Unevaluated. But now we meet another problem, xx will also be kept unevaluated, and thus Thread cannot be effective. So we have to put therein the evaluated List xx. We can do this by Substitution (Replace..), as an (evaluated) function argument or using With. I choose the latter: Evaluate[Fold[#2.#1 &, {0, 0, 1}, With[{xx = xx}, Thread[Unevaluated[evalmat[xx, #, z]]]]]] & /@ yy {{x3*z^5 + x3^2*(x1*x2^2*z + x2*y1^2*z + x2*z^3) + x3*y1*(y1^3*z + y1*z^3), y1*z^5 + y1^2*(y1^3*z + y1*z^3), z^6}, {x3*z^5 + x3^2*(x1*x2^2*z + x2*y2^2*z + x2*z^3) + x3*y2*(y2^3*z + y2*z^3), y2*z^5 + y2^2*(y2^3*z + y2*z^3), z^6}} Finally: if you meet conflicting goals (memory and speed), tackle the problem which hinders most (memory in your case). Best wishes -- Hartmut