Message Boards Message Boards

0
|
7601 Views
|
17 Replies
|
6 Total Likes
View groups...
Share
Share this post:

Solve a fraction when the denominator approaches zero?

Posted 7 years ago

I have a function (1+b)^2/Sqrt[1-b^2] that I need to solve as b->1, where b is real. The calculation exceeds the precision that Mathematica is using and I am not sure how to increase the precision if I use Solve...or even if this is the best way to get a solution. Is there another way of approaching the solution to this function from a purely analytical approach or using complex variables?

Regards,

Luther

POSTED BY: Luther Nayhm
17 Replies
Anonymous User
Anonymous User
Posted 7 years ago

I'll mention a lesson I experienced.

limited precision (ie on a calculator) can become complicated and need special handling even when arbitrary precision is used. on big problem with it is division.

an example is this: your program must use the result of a previous calculation to decide whether to go left or right on a calendar day / choice. if your precision is "perfect" you know theoretically it will go the right direction. HOWEVER: to divide mathematica must find factors (LCD/GCD etc, eventually), and that is where a flaw can enter

the answer is that some programs must be double careful to proceed in a manner where whether scanning forward or backward through digits never can land on ANY digit at ANY precision which causes their equation to work in one direction but not the other due to "extremely small error". this can be tricky. but your tolerance "at the turning point" cannot be theoretical upon division: you have to develop a way to track which side things will fall on (because ultimately: division is not a perfect thing, not way down deep)

in texas they offer courses on how to do scientific math with limited precision calculators (using non-junk calculators, TI calculators) for just such kinds of "programs" that approach a theory number that can't really be represented. one technique is insuring division isn't used. another is something like newton approx. where tolerance is known ahead and assured.

POSTED BY: Anonymous User

Ok, I might not understand the point of this discussion, but you appear not to be interested in

Limit[(1 + b)^2/Sqrt[1 - b^2], b -> 1, Direction -> "FromBelow"]

or

Limit[(1 + b)^2/Sqrt[1 - b^2], b -> 1, Direction -> "FromAbove"]

which is what Daniel suggests early on. Here is a different way of writing what Daniel has done in a much shorter way. First, I try to figure out for which "b" I get a certain value:

$MaxPrecision = 550    
sols = Solve[ff == c, b]

enter image description here

The first solution is the real valued solution. Now you can use:

sols[[1]] /. SetPrecision[c -> 9999999999999999999999999, 450]

and obtain something identical (plus/minus a couple of digits) as the simpler code before:

NSolve[ff - 9999999999999999999999999.`450 == 0, b]

You can get really close to the pole:

$MaxPrecision = 950
sols[[1]] /. SetPrecision[c -> 9999999999999999999999999999999999999999999999999999, 850]

enter image description here

I do not really have an easy means of verifying these numbers though. How close to the pole do you want to go?

There are not many applications, I am aware of, that need to such a high precision though, so I would like to hear of yours. Always nice to be able to tell students why one would need these high precisions "in real life".

Best wishes,

Marco

POSTED BY: Marco Thiel

This seems to be essentially the same what I noted above....... :)

POSTED BY: Hans Dolhaine

Oops, sorry, had not seen that.

BW,

Marco

PS: I'd still like to see the "application of this".

POSTED BY: Marco Thiel
Posted 7 years ago

I confess to using bad syntax and descriptive language, but you have essentially answered my questions, which I asked because I am only an intermittent user of the more advanced features...advanced to me...in Mathematica.

The application is an esoteric one of applying Einstein's radiation pressure model to a variety of scenarios. In the 1905 paper by Einstein, he describes how the Poynting Vector is impacted by detection in a moving frame...it is just a radar problem: what is the actual frequency that impacts a moving object? From Einstein's perspective, the radiation is Doppler shifted into the moving frame and then the interaction occurs.

I converted this model into a linear model using photons and flux and, ultimately, radiation pressure. Einstein's model omits what I have called a compression factor. This is just a geometrical factor that allows the number of photons that impact an object to depend on whether the object is moving toward a source or against the source of radiation. An example is that of a helium filled balloon held in a wind. When the balloon is stationary, the number of molecules of air hitting the balloon is constant. When I let the balloon go, it quickly accelerated to the wind speed and the number of wind molecules hitting the balloon drops nearly to zero. Likewise, if I walk into the wind holding the balloon, more molecules hit the balloon per unit time increasing the force on the balloon.

I generalized these factors for movement of an object obliquely through a radiation flux, but for many scenarios, we are moving directly into or against the radiation. I then began a parametric study of the effective radiation pressure on various sized objects moving at various speed in various directions. The question I posed was how to solve one of these scenarios, which is somewhat of a science fiction type, though I have also looked at solar sails, radiation sails in general, radiation forces on spacecraft, and particle accelerators. Some of the results are quite surprising and too quickly dismissed as having no practical consequence.

Thanks you for your patience and forbearance.

Luther.

POSTED BY: Luther Nayhm
Posted 7 years ago

Thank you. This is another approach I had not considered. The issue was about forcing the precision to be high enough that Mathematica actually gave me an answer. My post above indicates I was looking for the radiation pressure on a large object moving at light speed. The issue was, when does the up-shifted radiation cause sufficient radiation pressure to force the object to come to a terminal speed...assuming some constant thrust or driving force. The answer is that for a massive object, depending on the driving force, the object must reach ~25 nines the speed of light to reach a terminal velocity...but in the spirit of full disclosure, the object would be fried to a crisp at much lower speeds, perhaps below 0.9c.

Regards,

Luther

POSTED BY: Luther Nayhm

Welcome to Wolfram Community! Please make sure you know the rules: https://wolfr.am/READ-1ST

Please post your full code of what you are trying to do.

POSTED BY: EDITORIAL BOARD
Posted 7 years ago
ff = ((1 + b)^2)/Sqrt[1 - b*b]

(1 + b)^2/Sqrt[1 - b^2]

LogLogPlot[ff, {b, .9999999999999999999, .99999999999999999999}]

The limits of the plot are arbitrary here but need to be many nines larger for the circumstances I am looking for.

I originally tried Solve[ff-.99999999999999==0,b] and the solution was undefined by Mathematica. My solution should be at about 25 nines, but that is a guess at this point, since I cannot find an exact solution. If there is a non-Mathematica process that could be used, I am not familiar with it.

Luther

POSTED BY: Luther Nayhm

Please make sure you know read all rules: https://wolfr.am/READ-1ST

The rules explain how to format your code properly. If you do not format code, it may become corrupted and useless to other members. Please EDIT your posts and make sure code blocks start on a new paragraph and look framed and colored like this.

int = Integrate[1/(x^3 - 1), x];
Map[Framed, int, Infinity]

enter image description here

POSTED BY: EDITORIAL BOARD
Posted 7 years ago

The text I input was copied right out of the Mathematica worksheet, so it is the code I was executing.

It was not my intention to have the code executed. I was looking for information on if and whether Mathematica could be used to solve such a problem that borders on the pathological. I also did not want to open the discussion onto why I am trying to do this particular calculation.

POSTED BY: Luther Nayhm

Have you tried taking a limit as b approaches 1 from the left side?

POSTED BY: Daniel Lichtblau
Posted 7 years ago

Yes, since b<1. It is confusing since 0.99 is smaller than 0.999.

POSTED BY: Luther Nayhm

I'm not sure what is confusing about .99 being less than .999, or what that might have to do with taking a limit (which is what I had asked). I will add that the entire thread has been remarkably murky in terms of stating what exactly is wanted and what exactly has been tried. Let me comment on some related issues that I believe are being raised.

(1) Precision in a plot range. The input shown in an earlier post does not produce a plot (this fact, and the error message, should have been noted in that post).

ff = ((1 + b)^2)/Sqrt[1 - b^2]
LogLogPlot[ff, {b, .9999999999999999999, .99999999999999999999}]

(* During evaluation of In[122]:= Plot::plld: Endpoints for b in {b,0.9999999999999999999,0.99999999999999999999} must have distinct machine-precision numerical values.

Out[122]= LogLogPlot[ff, {b, 0.9999999999999999999, 
  0.99999999999999999999}] *)

What to do about it? One remedy is to widen the range. You can let the right bound be 1 for example. This next will give a picture that might be informative.

LogLogPlot[ff, {b, .9999, 1}, PlotRange -> All]

enter image description here

(2) Solving runs into precision problems. This can be alleviated by explicitly using higher precision in the input. Notice the back-tic in the input below. It parses the input number to 20 digits precision.

Solve[ff - .99999999999999`20 == 0, b]

(* Out[111]= {{b -> -5.00000000000006875*10^-15}} *)

(3) More confusion. If the interest is in b approaching 1, that's not the same thing as ff being close to 1. So the use of Solve seems unhelpful for the stated problem (as best I can interpret it). I would guess what is wanted is for ff to be large and then find b. This requires considerable precision, roughly twice as many digits in for the number you want out. Some analysis of the expression might help to explain this (e.g. expand in Puisiux series centered at b=1). Anyway, here is an example. Notice the decimal point is at the right rather than left side of the value for ff.

NSolve[ff - 99999999999999999999.`40 == 0, b]

(* Out[139]= {{b -> 0.9999999999999999999999999999999999999992}} *)
POSTED BY: Daniel Lichtblau
Posted 7 years ago

Yes...we are getting where I hopped to get. I do not have the background or vocabulary in Mathematica to be specific on what I am trying to do. You introduced some syntax that I was unfamiliar with but does get to the correct approach. I was only plotting "in the region of a solution" because I could not solve the equation with enough precision to get the answer I wanted.

In #2 above, b -> a negative value. It should be positive. Otherwise, you have solved my problem and answered my question.

If you are curious as to why I am making this calculation, let me know and I will attach a short document describing what I am looking for. It is a physics problem associated with modern technologies, which I would have posted in the physics section if I had wanted a discussion on the physics and how to use Mathematica to solve some technology modeling issues.

Thank you, though, for taking the time to answer this posting.

POSTED BY: Luther Nayhm

What do you mean by "solve"? Are you looking for the limiting value? Fi so, it blows up (that is, goes to plus or minus infinity depending on the direction of approach).

POSTED BY: Daniel Lichtblau
Posted 7 years ago

Sorry. I want to know what the values are as b->1 with b<1 but arbitrarily close to unity. I have a value I am looking for and wish to find b to solve the function for the large target value. In fact, the target value is so large that Mathematica will not solve the function because the denominator becomes so small.

POSTED BY: Luther Nayhm

I understand you have a function f given by the following equation e1

In[1]:= e1 = f == (1 + b)^2/Sqrt[1 - b^2]

Out[1]= f == (1 + b)^2/Sqrt[1 - b^2]

I square both sides and obtain

e2 = #^2 & /@ e1

f^2 == (1 + b)^4/(1 - b^2)

Having in mind that ( 1 - b^2 ) = ( 1 - b ) ( 1 + b ) you may guess that you will arrive at an equation with b^3. FullSimplify does the job

In[3]:= e3 = e2 // FullSimplify

Out[3]= (1 + b)^3/(-1 + b) + f^2 == 0

Now use

e4 = Solve[e3, b]

and you will get three equations for b (two with imaginary parts) which give you b when you plug in f. With Mathematica's numerical power this should give you an answer to your question.

For example

Block[{$MaxExtraPrecision = 500}, N[e4[[1]] /. f -> 10^80, 6000]]
POSTED BY: Hans Dolhaine
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract