mmcdara

7891 Reputation

22 Badges

9 years, 63 days

MaplePrimes Activity


These are replies submitted by mmcdara

@Carl Love @dnaviaux

 

Thanks Carl for having corrected me. 
As you've already noticed, I sometimes use the wrong term, and in this cas I'm all the more shameful that I lnew it perfectly well  "importance factor" or, as you wrote "factor that contributes to the error".
That's what happens if you think in your vernacular when writing in an other language.

In any case, here is @dnaviaux   a toy example which whows of to assess the effect of each factor.
Remember I use in a prcious reply the notional expression y=f(a, b, ...).
When all the (random variables which model) a, b, .. are mutually independent, the propagated variance writes
Variance(y) = diff(y, a)^2 * Variance(a) + diff(y, b)^2 * Variance(b) + ...

A classical definition of the "effect of a factor" (let's take "a" for instance) is 
Effect(a) = diff(y, a)^2 * Variance(a) / Variance(y) .

As SEA returns Variance(y), the method to obtain each "Effect" is simply to run SEA as many times as the number of factors while keeping fixed all the factors but the one of interest.
For instance:

  • ​​​​I write for short Variance(y) = SEA(a, b, c, ...)
  • Each a, b, c is defined this way : a := Quantity(ma, sa)), b := Quantity(mb, sb)), ...
  • Let a0 := Quantity(ma, 0)), b0 := Quantity(mb, 0)), ...
  • To assess the effect of factor a you just have to compute 
    SEA(a, b0, c0, ...) / SEA(a, b, c, ...)
     

SEA_SI.mw

 

PS 1 : I have said many times that the factors a, b, ... are modeled by random variables (when you use SEA for instance). This raise many questions:

  • If your Uncertainty is a constant value (lets say a = ma +/- Ua) then how are you going to define the value of sa when you define a as a:=Quantity(ma , sa)?
    • if you think the uncertainties are gaussian you will probably write something like 
      sa = Ua/2 or sa = Ua/3
    • if you believe they are uniform, then you will use sa = Ua/sqrt(3)
       
  • If your Uncertainty is a percentage (lets say a = ma *(1 +/- Ua)) this means the parameter  sin a:=Quantity(ma , sa) is not constant, which can be very complicated to handle if its mean (nominal) value ma varies in large proportions

PS 2 : If not all the funcertainties are mutually independent the formula for computing the effect of a factor is more complicated, but the more important point is that the definition of this effect is less simple.
The formula I gave above only accounts of what is called the "elementary effect" of (there) factor a.
If the uncertainties of a and b are correlated, then Variance(y) contains a term
2*diff(y, a)*diff(y, b)*Covariance(a, b).
This means that a contributes to Variance(y) through its interaction with b "interaction effect".
For a given factor there exist an "elementary effect" and a "total effect" wich includes the "elementary effect" plus all the effects through interactions this factor has with the other ones.

@ecterrab 

Thank you Edgardo,
Happy Holidays too

@dnaviaux 

Yes you can.
I have a worksheet I wrote one year ago for the needs of an internal training in my company.
I can provide it to you but it is presently in french.
What I propose is to translate it into English, or to build a toy example to show you can estimate the different contributions to the error on your quantities of interest.

The first step is probably to agree on the way to measure the effect of a contributor.
 

@Carl Love 

Thank you for lengthy Reply.

I'm still analyzing it, but to cut down on any polemic about the phrase "You already said me something like  "an answer is the definitive and unique answer to a question while a reply is just one among many proposals to answer this problem", I admit that I may have misinterpreted what you said or may have misattributed to you what others have said to me..
Please consider this a mistake on my part.

@nm 

In some sense it's better not to get a solution than to get a wrong one.

@dnaviaux 

Hi, 

What  ScientificErrorAnalysis  basically does is the following:

  • Let y=f(a, b, ...),
  • assume a, b, ... are random variables tthe only constraint is that the mean and the variance exist),
  • let X = < a, b, ...>, E(.) the operator "mathematical expectation" (the "mean" if you prefer) and let
    M = E(X) = < E(a), E(b), ...>,
  • a, b, .. are named "Quantity" in  ScientificErrorAnalysis,
  • then:
  1. E(y) ~ f(M) = f(E(a), E(b), ...))
  2. V(y) ~ E ( (X-M)^+ . H . (X-M)) 
    where X-M is the vector < a-E(a), b-E(b), ... >
    and H is the hessian matrix of f
  3. In the expression of V(y) appear tems of the form diff(y, a)^2*Varviance(a)
    and terms of the form 2*difff(y, a)*diff(y, b)*Covariance(a, b)


This is probably what you have done... even if the fact you mention the Jacobian of (here) y suggests that tou do something like Uncertainty(y) = diff(y, a)*Uncertainty(a) + diff(y, b)*Uncertainty(b) + ...

This is just a different point of view.
ScientificErrorAnalysis implements what is defined as "best practices" in the Guide to Uncertainty Measurements (GUM) gum.html 
See also here JCGM_100_2008_E.pdf, paragraph 5 for a detailed explanation of the "variance propagation formula" I talked about above.


Assessing the mean and variance of (here) y by these "propagation formulas" gives a good idea of their exact values if the function f can reasonnably be approximated by  a Taylor expansion of order 1 at point M
that is if f(X) ~ f(M) + grad(f(X))X=M . (X-M)+O((X-M)2).

But this is possible to derive higher order "propagation formulas", or evven "propagation formulas" for moments of higher orders (mean=moment of order 1, variance = (centered) moment of order 2).

My personal point of view, is that it's better to use ScientificErrorAnalysis than  Tolerance for the former is consistent withe the GUM's recomendations while the latter is not.
As a lot of people like tio have a final result of the form y=Y +/- E%, the GUM also gives recommendations to estimate the value of E:

  1. First use "propagation formulas" to assess Variance(y) (thus its standard deviation s(y))
  2. Secondly define a "reduced uncertainty" k (k usually taken to 2 or 3)
  3. Let's suppose s(y)=0.1 and k=2
    And say "E = 0.2 with a reduced uncertainty k=2" 

@vv 

Seems clear, thanks

@Carl Love @acer

"I don't know why you invariably submit Replies rather than Answers. There are buttons for each, just below a Question's body."

In April 19 2020 (Replace 'with' with 'uses') Carl wrote me this "Please stop posting things that should be considered Answers as Replies"

I understand the phrase "Please stop posting items that should be considered answers and not replies" but I was not sure to understand correctly Carl's phrase (in red, above).
This is why I made this reply the day after 

  • About Answers and Replies... could you explain clearly what difference you make between the twos (in french the both translate into the same word). You already said me something like  "an answer is the definitive and unique answer to a question while a reply is just one among many proposals to answer this problem". Is this right?
    With this in mind what I wanted to send was a reply, not an answer: if I sent an answer, please consider it as a mistake on my part and not as a deliberate act.
  • PS: copy-paste your phrase "Please stop posting things that should be considered Answers as Replies." into DeepL, translate it to french, take the result and translate it to english in order to understand the ambiguity between answer and reply.
    Maybe the term reply could be appropriately replaced by proposal?


I do not time to search the first blue phrase on Mapleprimes, perhaps I did not even repoduce it correctly?
But I always had it in mind when I dat to decide to "answer" or to "reply".
And because I do not have the pretension of saying the Truth (with capital T), or to give the most pertinent answer/reply, I always decided to play it modest and to send a non definitive contribution, so a reply.
And it is perfectly clear in the first item of my reply to Carl.

As I have never been contradicted on this point by anyone, I have made a point of answering only with "replies"

Now what you both say me is that the "default" attitude would be to send answers instead of replies?
There is definitely something I do not understand.
And I've just spent some minutes to read carefully the whole Mapleprimes' help without finding any explanation about the way to answer/reply a question.
The only hint is this text in the bottom left of the window "Be the first to answer this question." which suggests one must answer and not reply.

I thank you both for having changed "Replies to Answers, when appropriate"... but you see the ambiguity here: what does "appropriate" mean?
This suggests that, sometimes, a reply is the good choice... but who must decide of this, the one who answers/replies or a moderator?

@dnaviaux 

Thanks, I've loaded it.
It's rather complex and I need time to understant it clearly.
 

@jalal 

Default rendering in M&ple 2015 and 2020 are different.
Unfortunately I do not have this mlatter version at fome and I'm in hollidays for still 2 weeks, so I won't be able to help you.
Look to HighlitVertex for some clues.

Concerning notA for insttance, I guess you mean that the overbar is too close to A?
Try one of these possibilities (in the last one I use the conventional "negation" symbol)
 

notA := `#mover(mo(A),mo("&#95;")))`;
notA := `#mover(mo(A),mo("&#x305;")))`;
notA := `#mover(mo(A),munder(mo(""),mo("&#x305;")))`;
notA := `#mrow(mo("&#xac;"),mo(A))`;

Let me know

@Carl Love @dnaviaux

I agree with Carl.

In a few words:

  • Tolerances does basically what everybody has been taught, basically delta(a+b) = delta(a)+delta(b) (and a little mode complicated expression for functions of a and b)
  • ScientificErrorAnalysis does probabilistic uncertainty propagation : for the same toy example it would find Mean(a+b) = Mean(a)+Mean(b) and Variance(a+b)=Variance(a)+Variance(b)+2*Covariance(a, b)

All these relations come from first order or second order Taylor expansion of the "formula" (herea+b) you want to assess some uncertainty about.
This means none of these packages can give you an hint about the probability distribution of your quantity of interest (for instance the distribution of a+B given those of a and b... excepted if extremely simple cases).

Tolerances  is more an engineering tool to obtain a rough estimation of this uncertainty, while ScientificErrorAnalysis is more a tool devoted to statistical analysis.

Measurement devices or manufacturing tolerances are generally expressed as percentages of a reference value (for instance 10 ohms +/- 5%).
The first question is "how can I model this information?
If you use Tolerances  you will be able to define in a simple way these "+/- 5%", but if you use ScientificErrorAnalysis , that is if you adopt a statistical point of view, you need say for instance :

  • 5% represent 2 standard deviations of the (assumed) gaussian distribution of the manufacturing errors
    (why not 3 or 4?)
  • 5% represent the sqrt(3) standard deviations of the (assumed) uniform distribution of the manufacturing errors

Even if ScientificErrorAnalysis only requires means and variances to work, you see that the definition of the variances is not obvious.

None of these packages can handle common situations like

  • the measurement error is 0.1% in the range V1..V2
  • and 1% for V > V2

In this case (I assume exact computation cannot be done), only intensive simulation (Monte Carlo for instance) or ad hoc derivation of the propagated means and variances formulas, is reliable.

In the stationnary case of an electrical circuit I assume you have a linear system of the form Z*I=V where Z  is the impedance matrix, I a  vector of intensities and V a vector of potential differences?
This case will be "technically" handled either by Tolerances and ScientificErrorAnalysis  (with an ad hoc choice of the model of uncertainties).
The unsteady state ican be handled by  ScientificErrorAnalysis  (I've done it a lot of times but I cant' say if Tolerances would do the job).


About the generation of a report...
Here again I completely agree with Carl.
I use to generate LaTeX reports using printf but it is always in ad hoc way which depends on the report you want to write. So nothing that can be done completely automatically.
If your major concern is the writting of the report, perhaps Maple could not be the best choice:

  • I do not know Matlab enough, but Python (Jupiter NoteBooks) or R (Markdown documents) both enable writting a "LateX like" report in a rather simple way, with embedded code and/or results (even interactive graphs if you have the correct PDF version).
  • The equivalent of the package ScientificErrorAnalysis  exists in R, I do not know if the Matlab's Statistics Toolbox possesses such a feature?
  • Unfortunately there is no such facility in Maple ... (but I do not really know what can be done with DocumentTools or Embedded components)
     

Give us an example for us to see what you really want to do

Could you give an example of what you want to do?
 

@Kitonum 

The error was on the second line of Eqns := ...
Your line contained 

(xH(t)*cos(theta(t))+yH(t)*sin(theta(t)))/(x*H(t)^2+yH(t)^2)^2.5; 

instead of 

(xH(t)*cos(theta(t))+yH(t)*sin(theta(t)))/(xH(t)^2+yH(t)^2)^2.5; 

 

@Kitonum 

There is a function named H somewhere in your code that generates an error (probably a typo xH?, yH?)

PS:

  1. I use Maple 2015)
  2. Replacing Eqns by OP's Eqns no longer generates an error.


 

# before the call to dsolve
print(indets({Eqns, ICs}, function));

 /                      d             d             d          
{ H(t), cos(theta(t)), --- omega(t), --- theta(t), --- vxH(t), 
 \                      dt            dt            dt         

   d           d          d                             
  --- vyH(t), --- xH(t), --- yH(t), omega(0), omega(t), 
   dt          dt         dt                            

  sin(theta(t)), theta(0), theta(t), vxH(0), vxH(t), vyH(0), 

                                    \ 
  vyH(t), xH(0), xH(t), yH(0), yH(t) }
                                    / 

# error message
Warning, The use of global variables in numerical ODE problems is deprecated, and will be removed in a future release. Use the 'parameters' argument instead (see ?dsolve,numeric,parameters)
Error, (in DEtools/convertsys) the ODE system does not contain derivatives of the unknown function H

 

@Aragorn

In the loop "for n from 1 to periods do" (the syntax you use is unnecessarily complicated), write this line before trying to compute thetavalues.

print(plot('thetas(t)'-n*TH, t = 0 .. 10^6))

You will see that 'thetas(t)'-n*TH  is ptactically constant and eequal to -n*TH.
This is normal, n*TH is about 10^6 as thetas(t) is between 0 and 2*Pi by definition.
Thus you have no way to find a value of t such that thetas(t)=n*TH

 

Procedure orbit corrected (remove the red line as soon as you have fixed the problem I mentioned above)
orbit.mw

 

print(plot(['omegas(t)', 'thetas(t)', t = 0 .. 10^7], labels = ['omegas(t)', 'thetas(t)']))

First 93 94 95 96 97 98 99 Last Page 95 of 154