Dr. David Harrington

3481 Reputation

17 Badges

18 years, 93 days
University of Victoria
Professor or university staff
Victoria, British Columbia, Canada

Social Networks and Content at Maplesoft.com

Maple Application Center
I am a professor of chemistry at the University of Victoria, BC, Canada, where my research areas are electrochemistry and surface science. I have been a user of Maple since about 1990.

MaplePrimes Activity

These are replies submitted by dharr

@lcz Very nice. Looking through the paper, it seems that rather than take the pseudoinverse of L, you can take the regular inverse of L+(1/n)*J, where J is the matrix of 1's. (The wikipedia page seems to have missed this point.) Since you are not making use of GraphTheory's LaplacianMatrix, and making your own weighted version, you may as well generate the above matrix directly (called LJ in the code). Its diagonals are chosen so the row/column sums are 1. I'm assuming the regular inverse is more efficiently calculated than the Moore-Penrose pseudoinverse, though I didn't check. This leads to the code below (works for weighted and unweighted graphs).

  local W,M1,LJ,n,X;
  uses GraphTheory;
  if IsDirected(G) then error "Must be an undirected graph" end if;
  if IsWeighted(G) then W:=WeightMatrix(G) else W:=WeightMatrix(MakeWeighted(G)) end if;
  M1:=Matrix(n,n, (i,j)-> if W[i,j]<>0 then 1/n-1/W[i,j] else 1/n end if,'shape'='symmetric');
end proc:


@mmcdara @Carl Love Thanks for the appreciation. It seems like a deficit (bug?) that convert(..,LegendreP) here works while convert(..,LegendreQ) doesn't, though maybe some appropriate assumptions are needed. Since there were already complex arguments it wasn't obvious what they might be. I chose the change of variables so that z was the main argument of LegendreP. Ironically, if it hadn't half worked, I wouldn't have figured it out. I suppose that change of variables made the ode look more like the standard Legendre one.

Converting back does give slighly different hypergeoms, which no doubt are equivalent through some known relationships, but I'm not sure much can be concluded from this. 

You have edited your original question to be a different one, so the old answer no longer makes sense. To ask a related question in the same thread just reply to yourself. If it is a completely new question, start a new question. 

To have multiple plots of the same ode with different initial conditions use DEplot in the DEtools package. But you say you want a contour plot - you will need to specify more exactly what you want to get help with that.

@mskalsi Since you assigned to xi[1] and xi[2], xi is a table. So you cannot now use xi as a simple name. If you really want to do that, you need to use xi__1 and xi__2 (double underscore - also looks like a subscript). As @Axel Vogt says, eval works better than subs. But here to set (x^2+y^2)/2 to xi, you might have to use algsubs.

In terms of the specific question about how to plot a tetrahedron, you can supply a list of 4 points to tetrahedron in the plottools package. But that gets you a solid figure, and I think you want the lines. So the attached shows how to get lines beteeen all pairs of points and plot them. (You may have all the pairs already in your code, since those are the distances you wanted to be integer).


@tomleslie The simpler


works here. Of course I agree with your advice on using eval instead of assign.

@vv Thanks for the clarification.

@vv Nice. Perhaps I'm missing something, but it seems to me that this can be set up as a group with the set X as identity, and let Maple do the work. But it gives me an error message in my version (bug?), so I can't test it.

g:=CustomGroup(C,inverse=(a->X minus a),multiply=`union`);


@Tokoro n=10 plots in reasonable time and already is very dense so you can't see it well. So I'm not sure why you would want to plot the larger one. If you just want to do further calculations, you can just use the list of pairs of points in "lines", which is calculated rapidly.

I'm guessing it is the rendering of the plot which takes time - there are 26460 lines to plot. If so I don't think it can be made faster, since that part is outside the programmers control.

Edit: There is some inefficiency in the code, because I didn't know which plot you wanted. If you want only the lines, then just use 


but the time is still long because of the rendering.

@mehdibgh I don't know anything about this method, so can't comment further.

@mehdibgh I didn't attempt to follow the method, just fixed the Maple.

Edit: Q2 must surely be square for the congruences to work correctly.

@acer Thanks. The tubeplot was only intended as a workaround. I noticed the tubeplot flattening; it is the distortion because the plot is not constrained, and does indicate the inadequacy of the method.

You are right of course about the T,X,V,PX,PV:=rhs~(sol)[ ] ) line. I added the comment about ordering as a reminder, but eval is better (I may have used that in the past, but I don't remember.)

I was frustrated by the inefficiency of

fR := t->cols[signum(X(t)),signum(V(t))][1]:
fG := t->cols[signum(X(t)),signum(V(t))][2]:
fB := t->cols[signum(X(t)),signum(V(t))][3]:

and tried just 

f := t->cols[signum(X(t)),signum(V(t))]: with
color=f(t) or color=[f(t)[]]

but only got red (full evaluation?) and an error respectively.

I then debated producing the single number for the HUE instead, but I do have some time constraints.

I certainly hope the color function can be implemented into spacecurve; it seems a strange omission.

@mehdibgh M__11 should have been Mm__11 when creating RS. I always keep things real and symmetric, so use transpose instead of Hermitian transpose, and as I tried to emphasize before: tell the routines when things are symmetric (and positive definite). Now Eigenvectors works, is much faster and gives real eigenvalues. They still aren't accurate (according to your expected values), bocause your matrices are ill-conditioned. You could scale K so the values are more comparable to those in M, and then scale the resulting eigenvalues to compensate, but I don't think that will help much.


@mmcdara If you check the subscribe box when you post/answer/reply, then you will get email notifications and the thread will remain alive. But perhaps you won't notice this and find out.  :)

In terms of the possibilty of spacecurve coloring, when I found spacecurve gave me an error I assumed it couldn't be done. I was too lazy to check the plots,structure help page, but I used the word "think" not to be definitive.

So now I "think" only implicitplot3d's ISOSURFACEs can't take a color function. That one has frustrated me for a long time since I want to add phase colouring to those surfaces.

@mehdibgh Your Eig method doesn't work because the pseudoinverse is singular. 

I attach my version of your LV method, which avoids taking an inverse and makes some use of the fact that M is diagonal, though that could be avoided. I only did operations that kept the matrices symmetric, and avoided some of the slow conversions you had. Accuracy probably cannot be improved, because some data in M11 are small and inaccurate, so the inverse square root now has large inaccurate values. As you tried, some scaling can help but the matrices are all ill-conditioned. I'm guessing the input data does not have enough sig. figs.


5 6 7 8 9 10 11 Last Page 7 of 32