How To Solve Numerical

[Vella Massart]

Thisgives the correct answer in general, but appears to be more fragile than MATLAB. In the case I'm looking at now, it gives a very large negative number for x'/P*x as though there was overflow and wraparound or somesuch.

EDIT: Doing some experimentation I see that Eigen also fails to invert this matrix, where MATLAB has no problem. The matrix condition number is 3e7, which is not bad. Eigen gives the same (wrong) answer with the linear solve and a simple x.transpose() * P.inverse() * x. What's going on?

I would like to solve the differential equation in the document. However, Mathcad can only calculate the values if I use a positive factor. In Python the values are calculated without any problems. What am I doing wrong? The problem exists in Mathcad 15, as well as in Prime 6.0.

The root function is implemented in Mathcad so that it returns the real value, not the main root (the one with the smallest positive phase). Actually in math a power with a rational exponent is defined only for non-negative bases.

Now I have another challenge. If I solve the previous equation together in a differential equation system, Mathcad says "Convergence to one solution not possible, too many integrator steps" No matter which method or which step size I choose, the problem remains. Do you have another hint for me?

Hello All, I am relatively new to MathCAD and am having an issue developing a set of symbolically solved equations for a system. I have assigned values to the variables however throughout my analysis I would like to solve or simplify equations ONLY symbolically for a particular variable and KEEPING all variables as symbols, NOT replacing them with their numerical values.

These are all great options and if I am not mistaken they all tell me that if I want to solve purely symbolic sets of equations I have to rethink how I set up my document. In other words, I will be going through assembling equations for a variety of components that interact with each other (free body diagrams). I am going to get some relatively complex algebraic equations that describe the combined behavior.

It seems to me that the best option may be (1), i.e. to avoid assigning all numerical values to parameters and instead create a set of symbolic-only equations? Then once I have the overall solution equations I add the numerical values at the end and solve numerically?

Option 1) sure is a good choice, but I would recommend nonetheless turning all the calculations into functions depending on the variables you want to set in later - even if you do not intend to vary these variables.

Doing so may be troublesome but if you don't, the only way to use the symbolic results you calculated at the top would be to copy and paste the symbolic results which is not desirable and error prone.

In the case of quite elaborate and complicated equations, however, you run the risk of reaching the limits of Prime's capabilities, since the symbolics in Prime is relatively weak, especially compared to Maple or Mathematica.

Overall all these responses have cleared up a lot of basic functionality for me. One interesting thing that may be interesting for all is that in speaking also with PTC tech support, and having a little more complex equations going on looking for solution for a variable appearing in many locations, I was told that an upgrade to version 8 was required as the solver engine is much more robust, so I will be moving to that soon.

I think there is a system default for solution variable order. But control the order manually is difficult. If you want to write the variables any order what you want, use the bold = is one of the tips.

Thanks Werner. I thought your second solution might solve an annoyance I had but I couldn't get it to work. I'm solving a statistics equation that (to my knowledge) can only be solved symbolically because the factorials are otherwise too big for a computer to handle. The annoyance is I must evaluate my program symbolically, which just means the result must be shown. Only a big deal when you are evaluating dozens of results and it takes up pages (only 7 shown below). Ultimately I just want to plot the resulting values (probabilities) against the range values. So anyway can this be done without the results taking up so much page space? P7 attached.

My first idea was to use the modifier "float" so the symbolic result won't look that ugly and space consuming. But to my surprise this modifier has no effect! I also tried "return 1.0*prob" to force a float result, but to no avail. Surprisingly even "return 5*prob" did not return the expected vector multiplied by 5. I have no clue whats wrong with this.

So my second idea is to hide the ugly result in a collapsed area. In case you still want to show the calculation but hide the large result, you could use a function with a dummy parameter in front of a collapsed area, put the assignment with symbolic evaluation in a collapsed area and show the numeric result afterwards.

The function can, but must not have the same name as the variable and of course the name of the dummy argument can be chosen at will. You also may use a,b,x,y as function arguments, which sure would make sense, too.

Great work as always. As per your 2nd idea - I'd had the calculation off the page to hide it. I also had a dummy code as you suggested, but had just used the bold = so it wasn't calculated. Glad to see we went down the same paths.

As for combin, I didn't know that existed. Oddly (imho) it isn't listed under the statistics functions in MC, even thought it is under statistics in the online help menu. And you'd think it would show up if you searched for 'binomial', seeing as how most of the world calls it the binomial coefficient.

This is my first question on here. I'm trying to numerically solve the Schrdinger equation for the Woods-Saxon Potential and find the energy eigenvalues and eigenfunctions but I am confused about how exactly this should be done.

I've seen some mention of "tridiagonal matrices" being generated somehow, but am not sure what the elements of that matrix would be or how that applies to the problem. Leandro M. mentioned that "the discretization defines a finite dimensional (matrix) eigenvalue problem". This seems like the correct road I should be going down, but I haven't been able to find anything that explicitly explains this process or how the matrix is constructed. If this is the correct procedure, how is such a matrix constructed?

Numerov's method as described on Wikipedia is not how you want to proceed here. To give you an idea of how to proceed, let's start with a simplified version of the method. What I'm going to do is to just naively discretize the differential operator, like so:

and you may proceed to diagonalize it in whatever way you choose. Note that we call $A$ a tridiagonal matrix, for obvious reasons. You might want to take care of the boundary conditions first -- you do that by setting $\psi_1 = \psi_n = 0$ before you construct the matrix, which corresponds to setting the first and last columns to zero. This will give you some spurious eigenfunctions with zero eigenvalue that you can just discard. If you have different boundary conditions, you're out of luck -- I don't know of a way that makes it work in general.

This seems to be the way you're looking for but of course it's not the only way to do this. Griffiths describes one called "wag the dog" that consists of the following: you guess an initial value for an energy and compute the wavefunction however you want (RK, for instance). Chances are it won't be an eigenvalue so the wavefunction will blow up at infinity. It might go to $+\infty$ for large $x$ or it might go to $-\infty$. You now slowly vary the energy until the behavior at infinity "flips", that is, until the tail "wags". That will allow you to constrain the value of a single energy eigenvalue and give you the form of the wavefunction out to some maximum size that increases as the chosen $E$ approaches the exact eigenvalue.

Priliminaries:To get the eigenvalues from Numerov method you will need to know the wavefunction at the boundaries. Generally this would mean that you need to set the potential to infinity at the boundaries hence putting the wavefunction to zero at those points. For your potential, modify it as follows:

These kinds of problems easily can be solved by a new approach that I (Grant Bunker, Professor of Physics, IIT, [email protected]) have developed, that I call the Phase Method ( -method). It is fully described and explained in my new book. One of the many examples given is the Woods-Saxon potential. Working example code in Python and Mathematica is posted at that web site. The PM does not use discretization methods, which tend to give poor accuracy for discontinuous potentials, especially for states near the continuum.

The PM is more closely related to the Shooting Method (with the same level of accuracy) but does not require any guessing of eigenvalues or adjusting initial conditions. It automatically maps out the eigenvalues, and can provide wave functions and scattering phase shifts also. Instead of searching for convergent wavefunctions, as is the usual procedure, instead it finds bound state energy eigenvalues by rapidly counting the number of states below any specified energy (using something called the phase plot) and using that state counting function in an evolutionary-tree-like process to map out all the energy eigenvalues, preferably in parallel. It is particularly useful when there are nearly degenerate eigenvalues or you have no idea where they may lurk. Typical relative accuracy using extended precision arithmetic is 15-20 decimal digits. Once accurate eigenvalues are found, wavefunctions are easily calculated, perhaps surprisingly, without having to twiddle boundary conditions.

A simple serial Python implementation is given on the linked web site (as is full code in Mathematica) and is pasted below. The potential in this example is for linear v-shape (Abs[x]/2) but you can change it to anything you like, adjusting the range appropriately for the potential and upper energy limit that is selected. The x-range should extend well beyond the classical turning points at that energy, i.e. the points at which the total energy equals the potential energy. Full explanations for the unconventional aspects of this method are given in the book.Have fun!

3a8082e126