mm-392 === Subject: : multivariate optimization problemHow do I classify the following optimization problem:x - complex-valued vectorA_i - i-th matrixy_i - i-th objectiveMinimize in the least square sense a set ofquadratic formsx A_1 x* = y_1x A_i x* = y_2...x A_n x* = y_nwherex - complex-valued vectorA_i - i-th matrixy_i - i-th objectiveThank you=== === Subject: : Re: multivariate optimization problem> How do I classify the following optimization problem:> x - complex-valued vector> A_i - i-th matrix> y_i - i-th objective> Minimize in the least square sense a set of> quadratic forms> x A_1 x* = y_1> x A_i x* = y_2> ...> x A_n x* = y_n> where> x - complex-valued vector> A_i - i-th matrix> y_i - i-th objectiveThis is a nonlinear least squares problem.Introduce x=y+iz to reqrite it in real form, and you can pose itto a standard solver (Try NEOS).Arnold Neumaier=== === Subject: : Re: Confidence Intervals for the Mean of a Poisson Distribution by support1.mathforum.org (8.11.6/8.11.6/The Math Forum, $Revision: 1.9 primary) id i2I2E2R32018;>All,>I have found several references which quote the following for the C%>confidence interval of the mean of a poisson distibution from asingle>obeservation x.>( qchisq((100-C)/200, 2*x)/2, qchisq(1-(100-C)/200, 2*(x+1))/2 ) >Where qchisq(p,n) is the inverse of the one-tailed probability of the>chi-squared distribution with n degrees of freedom.>This result has been stated in several places without proof or>reference to a proof. I would like to find a proof or derivation of>this rather elegant result. Does anyone have any references or an>explantion?>--crazzellLook in Rohatgi VK's book: An introduction to probability theory andmathematical statistics. New York: John Wiley and Sons. Pages 212 to213.=== === Subject: : Savitzky-Golay MethodI have studied some paper on this topic, and the doubt that I have is howcan apply this method on a vector that has a spectrum signal?In other words, how to apply smoothing to the vector of data.=== === Subject: : Why mcc compiled windows exe is 20 times slower than m-files?? I have a matlab code that does a lot of matrix operations, thingslike A.*B.*C etc. in the main loop. The Windows EXE compiled with mcc runs 20 times SLOWER than them-file script. I used default options for mcc and mbuild, whichaccording to ther manual, has all optimizations turned on. Is this specific to the mcc compiler? Is there something I can do tospeed it up?=== === Subject: : Re: Numerical integration of Laguerre polynomials> what if you transform the integrals from [0,infy[ to ]0,1] > and use high order gauss-kronrod like quadpack does?Peter,thanks. I don't know much about Gauss-Kronrod so I'll have a look intothis.Martin Martin Stoll Institut f.9fr Theoretische Physik, Uni G.9attingen f: +49-551-39-7687 Tammannstrae 1, 37077 G.9attingen, Germany=== === Subject: : Re: Numerical integration of Laguerre polynomialsThe Laguerre polynomials L_{k}^{1/2} are in fact Hermite polynomials.(but perhaps you know this, because the problem comes from theharmonic oscillator area).Anyhow, we havesqrt(x) L_{k}^{1/2}(x) = (-1)^k 2^{-2k-1}/k! H_{2k+1}(sqrt{x}).So, if f and g allow, you can bring the u,v integrals tothe real line, with Gaussian weight functions, and youcan use Gaussian qudrature for Hermite polynomials(this is in fact a simpler problem, but will not helpyou further, I guess).The strange part is the x-integral over [-1,1].Can you give more details on f and g?============================================================ Nico M. Temme, http://homepages.cwi.nl/~nicot/C W I: Centrum voor Wiskunde en InformaticaKruislaan 413, NL-1098 SJ Amsterdam Tel +31 20 592 4240P.O. Box 94079, NL-1090 GB Amsterdam Fax +31 20 592 4199========================================================== =====> === Subject: : Numerical integration of Laguerre polynomials> for a while I've been stuck with the following problem and I was> hoping that someone here could give me some advise. It looks as if I'm> hitting the limits of double precision accuracy.> I need to do this triple integral numerically (physical background:> matrix element between harmonic oscillator states):> int_0^infty du sqrt(u) exp(-u) L_{kappa}^{1/2}(u) times> int_0^infty dv sqrt(v) exp(-v) L_{k}^{1/2}(v) times> int_{-1}^{+1} dx L_{kappa'}^{1/2}(f(u,v,x)) L_{k'}^{1/2}(g(u,v,x))> where L_{k}^{1/2} is an associated Laguerre polynomial. f and g are> simple functions. They are linear in x but not in u and v.> I use Gauss-Laguerre quadrature for the u- and v-integrals, and> Gauss-Legendre for x. This works very well for small indeces kappa,> kappa', k, k'.> However, already for kappa=0, kappa'=0 and k=30, k'=30, say, the> result of the triple quadrature becomes unreliable since the product> of the Laguerre polynomials in the integrand becomes very large> (1.D+300) whereas the quadrature weights beome very small (1.D-300)> and over/underflow occurs.> Neither could I write down an asymptotic formula for the integral> above.> Is there a better suited way to carry out the integral (I need> to it for many indeces...) in reasonable computational time?> Martin Stoll=== === Subject: : Re: Numerical integration of Laguerre polynomials> The Laguerre polynomials L_{k}^{1/2} are in fact Hermite polynomials.> (but perhaps you know this, because the problem comes from the> harmonic oscillator area).Nico,thank you for your suggestions. I was aware of the relation to theHermite polynomials (found it in the Abramowitz&Stegun) but I can'tsee if this could simplify the problem. Is there an advantage ofGauss-Hermite over Gauss-Laguerre quadrature?The x-integral over [-1,1] arises from an integration over a solidangle. Originally, it is the usual int_0^pi sintheta dtheta, thusx=costheta, and the reason for it's appearance in the Laguerrepolynomials is that their arguments are actually scalar products oftwo vectors.The functions f and g in my earlier post aref(u,v,x)=1/4 u + 3/4 v - 3/4 sqrt(uv) xg(u,v,x)=3/4 u + 1/4 v + 3/4 sqrt(uv) xso they're linear in x but not in u, v. Martin Martin Stoll Institut f.9fr Theoretische Physik, Uni G.9attingen f: +49-551-39-7687 Tammannstrae 1, 37077 G.9attingen, Germany=== === Subject: : Re: Numerical integration of Laguerre polynomials> for a while I've been stuck with the following problem and I was> hoping that someone here could give me some advise. It looks as if I'm> hitting the limits of double precision accuracy.> I need to do this triple integral numerically (physical background:> matrix element between harmonic oscillator states):> int_0^infty du sqrt(u) exp(-u) L_{kappa}^{1/2}(u) times> int_0^infty dv sqrt(v) exp(-v) L_{k}^{1/2}(v) times> int_{-1}^{+1} dx L_{kappa'}^{1/2}(f(u,v,x)) L_{k'}^{1/2}(g(u,v,x))> where L_{k}^{1/2} is an associated Laguerre polynomial. f and g are> simple functions. They are linear in x but not in u and v.Not sure if it will help but, writing f(u,v,x) = a + b x g(u,v,x) = c + d xyou could use a CAS to compute (and store) the x-integral in closed form for arbitrary a, b, c, d, and fixed kappa', k'. > I use Gauss-Laguerre quadrature for the u- and v-integrals, and> Gauss-Legendre for x. This works very well for small indeces kappa,> kappa', k, k'.> However, already for kappa=0, kappa'=0 and k=30, k'=30, say, the> result of the triple quadrature becomes unreliable since the product> of the Laguerre polynomials in the integrand becomes very large> (1.D+300) whereas the quadrature weights beome very small (1.D-300)> and over/underflow occurs.> Neither could I write down an asymptotic formula for the integral> above. > Is there a better suited way to carry out the integral (I need> to it for many indeces...) in reasonable computational time?Alternatively, would using the generating function LaguerreL[n, a, z] == SeriesTerm[(1 - t)^(-a - 1) Exp[(t z)/(t - 1)], {t, 0, n}]help. You can certainly obtain a simple closed form for the x-integral: Integrate[(E^((t (a + b x))/(t - 1))/(1 - t)^(3/2))* (E^((s (c + d x))/(s - 1))/(1 - s)^(3/2)), {x, -1, 1}] (E^(((c - d)*s)/(s - 1))*(-E^(((a - b)*t)/(t - 1)) + E^((2*d*s)/ (s - 1) + ((a + b)*t)/(t - 1))))/(Sqrt[1 - s]*Sqrt[1 - t]*(d*s*(t - 1) + b*(s - 1)*t))Cheers,PaulPaul Abbott Phone: +61 8 9380 2734School of Physics, M013 Fax: +61 8 9380 1014The University of Western Australia (CRICOS Provider No 00126G) 35 Stirling HighwayCrawley WA 6009 mailto:paul@physics.uwa.edu.au AUSTRALIA http://physics.uwa.edu.au/~paul=== === Subject: : Re: Numerical integration of Laguerre polynomialsPaul,thanks for your answers. I will have to think about the CASsuggestion. On the other hand, since Gauss-Legendre quadrature isexact for the x integral I'm not sure whether a termwise integrationand, afterwards, re-summation wouldn't lead to the same kind ofnumerical problem, that is, overflow?> Alternatively, would using the generating function> LaguerreL[n, a, z] == > SeriesTerm[(1 - t)^(-a - 1) Exp[(t z)/(t - 1)], {t, 0, n}]> help. You can certainly obtain a simple closed form for the x-integral:> Integrate[(E^((t (a + b x))/(t - 1))/(1 - t)^(3/2))*> (E^((s (c + d x))/(s - 1))/(1 - s)^(3/2)), {x, -1, 1}]> (E^(((c - d)*s)/(s - 1))*(-E^(((a - b)*t)/(t - 1)) + E^((2*d*s)/> (s - 1) + ((a + b)*t)/(t - 1))))/(Sqrt[1 - s]*Sqrt[1 - t]*(d*s*(t - 1)> + b*(s - 1)*t))I've also considered this possibility a while ago. The reason I didn'tproceed further was that I didn't know how to conveniently evaluatethe derivatives of this formula with respect to kappa' and k' such asto obtain the coefficient of the power series which gives me theoriginal integral I'm looking for. Or am I missing your point?Martin Martin Stoll Institut f.9fr Theoretische Physik, Uni G.9attingen f: +49-551-39-7687 Tammannstrae 1, 37077 G.9attingen, Germany=== === Subject: : Tributary Area AlgorithmI am trying to write a program to compute an area used for stressanalysis. The problem is as follows:1) I have a closed polyline that define an irregular shaped area in2D. This area mechanically defines a solid (YYY in picture below)2) Surrounding this area are 8 other irregular 2D closed polylinesthat define surrounding solid areas (xxx)3) Undefined areas around the 2D solid polyline areas (|||)picture: x = area defined by 2D polylines, | = areas not enclosed by2D polylinesxxx||xxx||xxx|||||||||||||xxx||YYY||xxx||||||||||||| xxx||xxx||xxx4) I need to find the closed polyline that surrounds the YYY area andis half-way between the YYY area and each of the xxx defined areas. The area defined by the half-way polyline defines a stress area thatthe YYY solid area must support.Any suggestions on the how to do this generically for the case ofirregular shaped solid YYY and xxx polylines?Rock=== === Subject: : Re: Tributary Area Algorithm> I am trying to write a program to compute an area used for stress> analysis. The problem is as follows:> 1) I have a closed polyline that define an irregular shaped area in> 2D. This area mechanically defines a solid (YYY in picture below)> 2) Surrounding this area are 8 other irregular 2D closed polylines> that define surrounding solid areas (xxx)> 3) Undefined areas around the 2D solid polyline areas (|||)> picture: x = area defined by 2D polylines, | = areas not enclosed by> 2D polylines> xxx||xxx||xxx> |||||||||||||> xxx||YYY||xxx> |||||||||||||> xxx||xxx||xxx> 4) I need to find the closed polyline that surrounds the YYY area and> is half-way between the YYY area and each of the xxx defined areas. > The area defined by the half-way polyline defines a stress area that> the YYY solid area must support.> Any suggestions on the how to do this generically for the case of> irregular shaped solid YYY and xxx polylines?> RockThe layout should be like this... YYY sould be in the center with thexxx on all sides.> xxx||xxx||xxx> |||||||||||||> xxx||YYY||xxx> |||||||||||||> xxx||xxx||xxx=== === Subject: : Re: Tributary Area Algorithm[Note F'up2 reduction to one group --- should have been done by OP!]> 4) I need to find the closed polyline that surrounds the YYY area and> is half-way between the YYY area and each of the xxx defined areas. This definition of what you're looking for appears too vague to beimplementable. A polyline or other curve can't be half-way betweentwo other curves (or polylines) just like that --- that position iswell-defined only among points or straight lines.You'll have to specify that more precisely.Hans- Broeker (broeker@physik.rwth-aachen.de)Even if all the snow were burnt, ashes would remain.=== === Subject: : Re: Tributary Area Algorithm> [Note F'up2 reduction to one group --- should have been done by OP!]Not sure what this first comment means...> 4) I need to find the closed polyline that surrounds the YYY area and> is half-way between the YYY area and each of the xxx defined areas. > This definition of what you're looking for appears too vague to be> implementable. A polyline or other curve can't be half-way between> two other curves (or polylines) just like that --- that position is> well-defined only among points or straight lines.> You'll have to specify that more precisely.Is it possible to include an image to show the geometry???Assume that you have an underground room with open spaces and solidpillars (much like a checker board) and assume that there is one solidpillar in the center, with 8 other solid pillars around the centerpillar. All around the solid pillars is a void space, i.e. air. Iwould like to find a polyline (a group of connected points) thatdefine an area around the center pillar. The points used to definethe closed polyline are obtained by determining the closest distancebetween the edge of the center pillar and the edge of an outer pillarand then locating the point half way between the two pillars.The real world application of this is to estimate the stress acting onthe center pillar by assuming the pillar must support an area that isenclosed by a polyline drawn around the center pillar with pointsmid-way between the center pillar and the outer pillars.I hope this clarifies the problem...It seems to me that one way to do this problem is to step around thepolyline casting rays to find the closest point. There are problemswith irregular pillars and when defining distances for outer pillarsthat are diagonal to the center pillar.Are there any other bright ideas for working this problem?=== === Subject: : Re: Ability In MathematicsX-NewsReader: GRn 3.2n February 9, 1999> So,given all this, have you any ideas as to why,for some people,> Mathematics ability comes so naturally whilst other people find it so> difficult?Clearly some people see patterns in Mathematics whilst> others,NO MATTER HOW MUCH THEY TRY,cannot.> Short answer -- different learning styles.And another short answer is: attitude.> I saw two examples of this when I was in university {engineering student}.[snip of very good stories]As an engineer, I've always been interested in math asa tool to accomplish a task. The mathematicians seemath as an end in itself. I couldn't care less aboutwhether a math process/concept can be proven. Does itprovide useful results? Therefore, I am not a goodmathematician in the sense of pure math, but I canappreciate the basis for the math that I use.Another short answer is the person presenting the material.I just barely made it through a diffeq course taught by agrad student reputed to be brilliant mathematician. Hebreezed through the material which was trivial to him, buthe boggled my mind to the point that I didn't evenrecognize a garden-variety first order equation on a testwhen I was using the same equations every day in my EE courses.=== === Subject: : Re: Keeping quaternion normalized> what's the recommended method of keeping the unit quaternion> representing orientation in a rigid body model normalized if the state> equation is integrated by some numerical method?> Zipfel, Modeling and Simulation of Aerospace Vehicle Dynamics,> recommends adding a term > k ( 1 - || q || ) q> to the RHS of the state equation, which effectifly uses the ODE solver> to move q back to the surface of the unit sphere as soon as it's> deviating. But what are the reasons against the simple method of just> setting> q <- 1 / sqrt( || q || ) q> after each integration time step?If you use a penalty term, you can use any high-performanceautomatic integrator from the shelf, while if you rescale aftereach step you need 1. to write your own code, and 2. are limitedto single-step (Runge-Kutta type) methods.Arnold Neumaier=== === Subject: : Re: Keeping quaternion normalized> You're thinging q <- q/||q||, right?No, actually I meant the sqrt. By ||q|| I mean the Cayley norm, whichis q_0^2 + q_1^2 + q_2^2 + q_3^2, contrary to the usual notions forEuclidian vectors.> Adding a term k(1 - sum(q_i ^ 2))q has the advantage that you don't have to> take a square root -- doing the four multiplies is bad, but taking a square> root is really going to add to your execution time.No it won't. I my application, the state integration takes up only anegligible part of CPU. My question aims more at possible hiddennumerical issues!Kind regards-GerhardGerhard Wesp o o Tel.: +41 (0) 43 5347636Bachtobelstrasse 56 | http://www.cosy.sbg.ac.at/~gwesp/CH-8045 Zuerich _/=== === Subject: : Re: Keeping quaternion normalized> You're thinging q <- q/||q||, right?> No, actually I meant the sqrt. By ||q|| I mean the Cayley norm, which> is q_0^2 + q_1^2 + q_2^2 + q_3^2, contrary to the usual notions for> Euclidian vectors.> Adding a term k(1 - sum(q_i ^ 2))q has the advantage that you don't haveto> take a square root -- doing the four multiplies is bad, but taking asquare> root is really going to add to your execution time.> No it won't. I my application, the state integration takes up only a> negligible part of CPU. My question aims more at possible hidden> numerical issues!> Kind regards> -Gerhard> -- > Gerhard Wesp o o Tel.: +41 (0) 43 5347636> Bachtobelstrasse 56 | http://www.cosy.sbg.ac.at/~gwesp/> CH-8045 Zuerich _/Shot down on all counts! Oh well.I strongly suspect that the original reason for servoing the length _was_calculation time, because the sqrt and divide operation can both be veryexpensive compared to multiplies and adds. I'd verify that adding the sqrtoperator (and the divide, for that matter) didn't slow things down.=== === Subject: : Re: Keeping quaternion normalized> Shot down on all counts! Oh well.Sorry, wasn't my intention :-)May well be that originally computation time was a problem, but nowadays this should no longer be the case.-GerhardGerhard Wesp o o Tel.: +41 (0) 43 5347636Bachtobelstrasse 56 | http://www.cosy.sbg.ac.at/~gwesp/CH-8045 Zuerich _/=== === Subject: : Re: Keeping quaternion normalized> Shot down on all counts! Oh well.> Sorry, wasn't my intention :-)> May well be that originally computation time was a problem, but nowadays> this should no longer be the case.> -Gerhard> -- > Gerhard Wesp o o Tel.: +41 (0) 43 5347636> Bachtobelstrasse 56 | http://www.cosy.sbg.ac.at/~gwesp/> CH-8045 Zuerich _/I'm pretty close to doing some work with this sort of math on processorswithout floating point units, in which case it will still be the case forme! Desktop Pentiums may not have this issue, but as soon as you put apower restriction on the processor you find that the world shifts.=== === Subject: : Re: Help, if you can>What do these numbers have in common - 7, 3, 10, 11, and 12 ??> One of the four words in Help, if you can shares the property,> as does one of the seven words in What do these numbers have in common======================================================== ====================Possible answers :Solution 1.===Let a_1=-7 , a_2=3 , a_3=10 , a_4 =11 , a_5= 12 , and (x-2)(x-3)(x-4)(x-5)L_1(x)= ---------------------- 24 (x-1)(x-3)(x-4)(x-5)L_2(x)= - ---------------------- 6 (x-1)(x-2)(x-4)(x-5) L_3(x)= --------------------- 4 (x-1)(x-2)(x-3)(x-5)L_4(x)= - -------------------- 6 (x-1)(x-2)(x-3)(x-4)L_5(x)= ---------------------- . 24Let us consider the polynomial Curve (C) of explicit equation(#) y=P(x)= L_1(x)*(-7) +L_2(x)*3 +L_3(x)*10 +L_4(x)*11 +L_5(x)*12 , x in R .Then the points P_k from R^2 P_k= (k, a_k) , k in {1,2,3,4,5}are on the graphic of polynomial curve (fourth degree) y=P(x) , x in R . However this is not the unique caracterization of {a_1,a_2,a_3,a_4,a_5}. =====Solution (II). Other property satisfied by these numers is :consider an arbitrary system of positive numbers W={w_1,w_2,w_3,w_4,w_5) ,select a real number p,p=/=0, and define SUM_{k=1 to k=5}w_k*|L_k(x)|^{p} *a_kR(x)=R(W;p;x)= ------------------------------------- . SUM_{k=1 to k=5}w_k*|L_k(x)|^{p} It's easy to see that R(j)= a_j , for j in {1,2,3,4,5} .Therefore P_k , k in {1,2,3,4,5} are on the graphic of (rational curve) y= R(x) , x in R .Perhaps help, Alex. ======================== === Subject: : tangent singularitiesBecause of useage of the tan(gent) rountine near itssingularity pi/2 (and odd multiples) I wanted checkwhat would happen at M_PI_2 with the following code:#include )which I discovered worked without any thought to thevalue of arg. I first thought it might be the C99concept of infinite numbers which would work finein this case (unless was zero ??!!).All docs for the tangent routine make comment aboutpotential failure, but so far I can't make it fail.Thoughts.=== === Subject: : Averaging quaternionsDoes anyone know if it is possible to take an average of regularlysampled quaternions to get a mean orientation (i.e. a mean rotationmatrix)? I seem to recall there being a trick involved but beyondre-normalizing the resuling (averaged) quaternion, I cannot rememberwhat it is.Cheers,Graham=== === Subject: : Re: Averaging quaternions> Does anyone know if it is possible to take an average of regularly> sampled quaternions to get a mean orientation (i.e. a mean rotation> matrix)? I seem to recall there being a trick involved but beyond> re-normalizing the resuling (averaged) quaternion, I cannot remember> what it is.I am sure I will be scolded by somebody, but I believe that you canaverage the quaternion components, and then normalize as you say.This is assumes that you are noise dominated. Also, there is one trick that I can think of, which is thatquaternions are degenerate. For each unique rotation, there are twopossible quaternions whose components have opposite signs. This isbecause a positive rotation about axis V is identical to a negativerotation about axis -V.If your system is capable of both signs indiscriminately, then youmust make the sign conventions uniform. For example, by always makingone component positive.Some advertising since you crossposted on the IDL newsgroup: I do havea fairly comprehensive quaternion IDL library on my web page.CraigP.S. http://cow.physics.wisc.edu/~craigm/idl/idl.html (under Math)--------------------------------------------------------- -----------------Craig B. Markwardt, Ph.D. EMAIL: craigmnet@REMOVEcow.physics.wisc.eduAstrophysics, IDL, Finance, Derivatives | Remove net for better response------------------------------------------------------ --------------------=== === Subject: : About Fortran complie and buildI used Compag Visual Fortran to compile some soucefilem. But it wasinterupted by the error and gave the follow cue:*******--------------------Configuration: mod_common - Win32Debug--------------------Linking...dfor.lib(DFORMAIN.OBJ) : error LNK2001: unresolved external symbol_MAIN__Debug/mod_common.exe : fatal error LNK1120: 1 unresolved externalsError executing link.exe.mod_common.exe - 2 error(s), 0 warning(s)*******I donnt know why?