Calculer la moyenne moyenne de deux nombres

41

disclaimer: la moyenne est composée par moi

Définissez la moyenne arithmétique de nombres sous la forme Définissez la moyenne géométrique de nombres sous la forme Définissez la moyenne harmonique de nombres comme M _ {- 1} (x_1, ..., x_n) = \ frac {n} {\ frac {1 } {x_2} + \ frac {1} {x_2} + ... + \ frac {1} {x_n}} Définissez la moyenne quadratique de n nombres comme M_2 (x_1, ..., x_n) = \ root \ of {\ frac {x_1 ^ 2 + x_2 ^ 2 + ... + x_n ^ 2} {n}} La moyenne ( M_M ) est définie comme suit: Définissez quatre séquences ( a_k, b_k, c_k, d_kn

M1(X1,...,Xn)=X1+X2+...+Xnn
n
M0(X1,...,Xn)=X1X2...Xnn
n
M-1(X1,...,Xn)=n1X2+1X2+...+1Xn
n
M2(X1,...,Xn)=X12+X22+...+Xn2n
MMunek,bk,ck,rék ) comme
une0=M1(X1,...,Xn),b0=M0(X1,...,Xn),c0=M-1(X1,...,Xn),ré0=M2(X1,...,Xn),unek+1=M1(unek,bk,ck,rék),bk+1=M0(unek,bk,ck,rék),ck+1=M-1(unek,bk,ck,rék),rék+1=M2(unek,bk,ck,rék)
 Les quatre séquences convergent vers le même nombre,MM(X1,X2,...,Xn) .

Exemple

La moyenne moyenne de 1 et 2 est calculée comme suit: commencez par Alors Le calcul ultérieur des séquences doit être clair. On peut voir qu'ils convergent vers le même nombre, environ 1,45568889 .

a0=(1+2)/2=1.5,b0=1∗2=2≈1.4142,c0=211+12=43≈1.3333,d0=12+222=52≈1.5811.
une1=1,5+1.4142+1,3333+1,58114≈1,4571,b1=1,5*1.4142*1,3333*1,58114≈1,4542,c1=411,5+11.4142+11,3333+11,5811≈1,4512,ré1=1,52+1.41422+1,33332+1,581124≈1.4601.
1,45568889

Défi

Étant donné deux nombres réels positifs, et ( ), calculez leur moyenne moyenneunebune<bMM(une,b) .

Cas de test

1 1 => 1
1 2 => 1.45568889
100 200 => 145.568889
2.71 3.14 => 2.92103713
0.57 1.78 => 1.0848205
1.61 2.41 => 1.98965438
0.01 100 => 6.7483058

Remarques

  • Votre programme est valide si la différence entre sa sortie et la sortie correcte n’est pas supérieure à 1/100000 de la valeur absolue de la différence entre les nombres saisis.
  • La sortie devrait être un nombre unique.

C'est du code-golf , alors le code le plus court gagne!

code-golf math mon pronom est monicareinstate
la source
Lien Sandbox
mon pronom est monicareinstate
Quelque peu lié: Calculez la moyenne arithmétique – géométrique
user202729
11
À quel point sommes-nous censés être précis?
Incarnation de l'Ignorance
2
Étroitement liée
Giuseppe
1
Pouvons-nous supposer que la première entrée est toujours plus petite que la seconde, comme dans tous vos cas de test? (Si ce n'est pas le cas, je vais revenir en arrière sur ma réponse en Java.)
Kevin Cruijssen Le

Réponses:

14

Wolfram Language (Mathematica) , 52 octets

#//.x_:>N@{M@x,E^M@Log@x,1/M[1/x],M[x^2]^.5}&
M=Mean

Essayez-le en ligne!

Dans ma première approche, j’ai utilisé ces commandes intégrées
Mean GeometricMean HarmonicMeanetRootMeanSquare

Voici quelques substitutions pour économiser des octets

HarmonicMean-> 1/Mean[1/x] par @Robin Ryder (3 octets enregistrés)
GeometricMean-> E^Mean@Log@xpar @A. Rex (2 octets enregistrés)
RootMeanSquare-> Mean[x^2]^.5par @A. Rex (4 octets enregistrés)

Enfin, nous pouvons affecter Meanà M(comme proposé par @ovs) et économiser 5 octets supplémentaires.

J42161217
la source
Économisez 2 octets en recodant GeometricMean
Robin Ryder
@RobinRyder Je crois que tu veux dire Harmonic .. nice!
J42161217
1
Économisez 8 autres octets :#//.x_:>N@{Mean@x,E^Mean@Log@x,1/Mean[1/x],Mean[x^2]^.5}&
A. Rex
@ovs edited .....
J42161217 Le
10

R, 70 69 67 octets

x=scan();`?`=mean;while(x-?x)x=c((?x^2)^.5,?x,2^?log2(x),1/?1/x);?x

Essayez-le en ligne!

-1 octet avec un meilleur conditionnement.
-2 octets en basculant en base 2.

Comme d'autres réponses, celui - ci utilise l'expression de la moyenne géométrique comme une moyenne arithmétique sur l'échelle logarithmique (ici en base 2):

M0(X1,…,Xn)=2M1(bûche2⁡X1,…,bûche2⁡Xn).

Il utilise également le fait que ∀k,rék≥unek≥bk≥ck , c’est-à-dire rék=max(unek,bk,ck,rék) . La condition unek=bk=ck=rék est donc équivalente à rék=M1(unek,bk,ck,rék) , qui est ce queutilise dans la boucle tandis que; ceci est réalisé en abusant de la syntaxewhiledont ne prend en compte que le premier élément lorsque la condition est un vecteur, d'où l'ordre dans lequel les moyens sont stockés. (Notez que nous pourrions aussi utiliserck car c'est le minimum des quatre, mais nous ne pourrions pas utiliserunek ou unbk dans la condition.)

Lorsque nous sortons de la boucle while, xest un vecteur constant. La finale ?xcalcule sa moyenne pour le réduire à un scalaire.

Robin Ryder
la source
1
Ne devrait-il pas être au lieu de l o g x n ? lnXnlogXn
Tau
@Tau Oui, je désignais logarithme naturel par , valeur par défaut de R. Quoi qu'il en soit, je l'ai maintenant remplacé par le logarithme en base 2 pour -2 octets. log
Robin Ryder
6

J , 34 octets

(31 en tant qu'expression sans l'affectation de variable f)

f=:1{(^.z,%z,*:z,[z=:(+/%#)&.:)^:_

Essayez-le en ligne!

Pour les fonctions aet b, a &.: b("un sous b" ( défi lié )) équivaut à (b inv) a b- appliquer b, puis a, puis l'inverse de b. Dans ce cas, la moyenne géométrique / harmonique / quadratique est la moyenne arithmétique "sous" logarithme, inversion et carré, respectivement.

utilisateur202729
la source
5

TI-BASIC, 42 35 34 octets

-1 octet grâce à @SolomonUcko

While max(ΔList(Ans:{mean(Ans),√(mean(Ans²)),mean(Ans^-1)^-1,e^(mean(ln(Ans:End:Ans(1

L'entrée est une liste de deux entiers dans Ans.
La sortie est stockée dansAns et est automatiquement imprimée à la fin du programme.

Les formules utilisées pour les moyennes géométriques, harmoniques et quadratiques sont basées sur l'explication de user202729 .

Exemple:

{1,2
           {1 2}
prgmCDGFB
     1.455688891
{100,200
       {100 200}
prgmCDGFB
     145.5688891

Explication:
(des nouvelles lignes ont été ajoutées à des fins de clarification. Elles n'apparaissent pas dans le code.)

While max(ΔList(Ans           ;loop until all elements of the current list are equal
                              ; the maximum of the change in each element will be 0
{                             ;create a list containing...
 mean(Ans),                   ; the arithmetic mean
 √(mean(Ans²)),               ; the quadratic mean
 mean(Ans^-1)^-1,             ; the harmonic mean
 e^(mean(ln(Ans               ; and the geometric mean
End
Ans(1                         ;keep the first element in "Ans" and implicitly print it

Remarques:

TI-BASIC est un langage à jeton. Le nombre de caractères ne correspond pas au nombre d'octets.

e^(est ce jeton d'un octet.

^-1est utilisé pour ce jeton d'un octet.
J'ai opté pour l'écriture ^-1car le jeton ressemble ֿ¹à un bloc de code.

√(est ce jeton d'un octet.

ΔList(est ce jeton de deux octets.

Tau
la source
Je pense que vous pouvez économiser une parenthèse en mettant la moyenne géométrique en dernier.
Solomon Ucko le
@SolomonUcko ah, merci de l'avoir remarqué! Je n'ai pas pensé à ça avant.
Tau
max(DeltaList(Ans-> variance(Ans.
lirtosiast
5

Java 10, 234 229 214 211 215 206 203 196 196 180 177 octets

a->{for(;a[1]-a[0]>4e-9;){double l=a.length,A[]={0,0,0,1};for(var d:a){A[2]+=d/l;A[3]*=Math.pow(d,1/l);A[0]+=1/d;A[1]+=d*d;}A[0]=l/A[0];A[1]=Math.sqrt(A[1]/l);a=A;}return a[0];}

-5 octets grâce à @PeterCordes .
-15 octets supplémentaires grâce à @PeterCordes , inspiré de la réponse R de @RobinRyder .
+4 octets parce que j'ai supposé que les entrées sont pré-commandées.
-27 octets grâce à @ OlivierGrégoire .

Essayez-le en ligne.

Explication:

a->{                        // Method with double-array parameter and double return-type
  for(;a[1]-a[0]            //  Loop as long as the difference between the 2nd and 1st items
                >4e-9;){    //  is larger than 0.000000004:
    double l=a.length,      //   Set `l` to the amount of values in the array `a`
           A[]={0,0,0,1};   //   Create an array `A`, filled with the values [0,0,0,1]
    for(var d:a){           //   Inner loop over the values of `a`:
      A[2]+=d/l;            //    Calculate the sum divided by the length in the third spot
      A[3]*=Math.pow(d,1/l);//    The product of the power of 1/length in the fourth spot
      A[0]+=1/d;            //    The sum of 1/value in the first spot
      A[1]+=d*d;            //    And the sum of squares in the second spot
    }                       //   After the inner loop:
                            //   (the third spot of the array now holds the Arithmetic Mean)
                            //   (the fourth spot of the array now holds the Geometric Mean)
    A[0]=l/A[0];            //   Divide the length by the first spot
                            //   (this first spot of the array now holds the Harmonic Mean)
    A[1]=Math.sqrt(A[1]/l); //   Take the square of the second spot divided by the length
                            //   (this second spot of the array now holds the Quadratic Mean)
    a=A;                    //   And then replace input `a` with array `A`
  }                         //  After the outer loop when all values are approximately equal:
  return a[0];}             //  Return the value in the first spot as result
Kevin Cruijssen
la source
En C, vous pouvez f+=Math.abs(d-D)<1e-9;obtenir une conversion implicite d'un résultat de comparaison booléen en un entier 0/1, puis double. Java a-t-il une syntaxe compacte pour cela? Ou est-il possible de le faire f+=Math.abs(d-D)et de vérifier ensuite que la somme des différences absolues est suffisamment petite ?
Peter Cordes le
1
Yup, pour vos cas de test, f>1e-8fonctionne comme une condition de boucle: 229 octets. a->{for(double f=1,D,A[],l;f>1e-8;a=A){D=a[0];A=new double[]{f=0,1,0,0};for(var d:a){f+=Math.abs(d-D);A[0]+=d;A[1]*=d;A[2]+=1/d;A[3]+=d*d;}A[0]/=l=a.length;A[1]=Math.pow(A[1],1/l);A[2]=l/A[2];A[3]=Math.sqrt(A[3]/l);}return a[0];}. Avec 1e-9, il est plus lent (environ deux fois le temps de calcul), obligé de faire plus d’itérations pour obtenir essentiellement 4 * d-Daussi petit. Avec 1e-7, c'est à peu près la même vitesse que 1e-8. Avec 1e-6, certains des derniers chiffres diffèrent d’un cas à l’autre.
Peter Cordes
1
La réponse de @ RobinRyder indique que la moyenne quadratique est toujours la plus grande et harmonique toujours la plus petite. Vous pouvez donc peut-être abandonner fentièrement et vérifier uniquement a[3]-a[2]<4e-9.
Peter Cordes le
1
@PeterCordes l==2||vous voulez dire ( joué au golf l<3|). Mais oui, bon point; Je l'ai ajouté. :)
Kevin Cruijssen le
2
180 octets en agrégeant des réducteurs agrégables.
Olivier Grégoire le
3

Charbon de bois , 40 octets

Ｗ‹⌊θ⌈θ≔⟦∕ΣθＬθＸΠθ∕¹Ｌθ∕ＬθΣ∕¹θ₂∕ΣＸθ²Ｌθ⟧θＩ⊟θ

Essayez-le en ligne! Le lien est vers la version verbeuse du code. Prend les entrées sous forme de tableau de nombres. Explication:

Ｗ‹⌊θ⌈θ

Répétez l'opération tant que le tableau contient différentes valeurs ...

≔⟦....⟧θ

... remplace le tableau par une liste de valeurs:

∕ΣθＬθ

... la moyenne...

ＸΠθ∕¹Ｌθ

... la moyenne géométrique ...

∕ＬθΣ∕¹θ

... le moyen harmonique ...

₂∕ΣＸθ²Ｌθ

... et la racine signifie carré.

Ｉ⊟θ

Convertit un élément du tableau en chaîne et l'imprime implicitement.

Neil
la source
3

PowerShell , 182 180 183 octets

$f={$a=$c=$d=$n=0
$b=1
$m=[math]
$args|%{$n++
$a+=$_
$b*=$_
$c+=1/$_
$d+=+$_*$_}
$p=($a/$n),$m::pow($b,1/$n),($n/$c),$m::sqrt($d/$n)|%{$m::round($_,9)}
if($p-ne$p[0]){$p=&$f @p}$p[0]}

Essayez-le en ligne!

Andrei Odegov
la source
171 octets
mazzy
@mazzy, impressionnant :)
Andrei Odegov
3

05AB1E , 26 24 23 octets

Δ©ÅA®.²ÅAo®zÅAz®nÅAt)}н

Essayez-le en ligne ou voyez les étapes de tous les cas de test .

-1 octet grâce à @Grimy .

Alternative de 23 octets pour la moyenne géométrique:

Δ©P®gzm®ÅA®zÅAz®nÅAt)}н

Essayez-le en ligne ou consultez les étapes de tous les cas de test .

Explication:

Δ         # Loop until the list no longer changes:
 ©        #  Store the current list in variable `®` (without popping)
          #  (which is the implicit input-list in the first iteration)
          #  Arithmetic mean:
  ÅA      #   Builtin to calculate the arithmetic mean of the list
          #  Geometric mean:
  ®.²     #   Take the base-2 logarithm of each value in the list `®`
     ÅA   #   Get the arithmetic mean of that list
       o  #   And take 2 to the power of this mean
          #  Harmonic mean:
  ®z      #   Get 1/x for each value x in the list `®`
    ÅA    #   Get the arithmetic mean of that list
      z   #   And calculate 1/y for this mean y
          #  Quadratic mean:
  ®n      #   Take the square of each number x in the list from the register
    ÅA    #   Calculate the arithmetic mean of this list
      t   #   And take the square-root of that mean
  )       #  Wrap all four results into a list
}н        # After the list no longer changes: pop and push its first value
          # (which is output implicitly as result)
Kevin Cruijssen
la source
23:Δ©P®gzm®ÅA®zÅAz®nÅAt)}н
Grimmy le
@ Grimy Merci! Je ne peux pas croire que je n'avais pas pensé à utiliser la longueur au lieu de Y2/4. :)
Kevin Cruijssen le
1
Un autre 23 qui améliore montre la similitude de la moyenne géométrique pour les autres: Δ©ÅA®.²ÅAo®zÅAz®nÅAt)}н. Malheureusement, il ne semble pas que nous puissions refactoriser tous ces ÅAproblèmes.
Grimmy
@ Grimy Oh, j'aime bien cette deuxième version. :) EDIT: Oups .. merci d'avoir remarqué mon erreur dans l'explication ..>.>
Kevin Cruijssen
Je ne programme pas très bien dans 05ab1e, mais pouvez-vous calculer des sommes et les diviser ensuite par la longueur plus tard?
mon pronom est monicareinstate le
2

Gelée , 25 à 24 octets

Wẋ4¹ÆlÆeƭ²½ƭİ4ƭÆm$€⁺µÐLḢ

Essayez-le en ligne!

Explication

                    µÐL | Repeat until unchanged:
W                       |   Wrap as a list
 ẋ4                     |   Copy list 4 times
                   ⁺    |   Do twice:
                 $€     |     For each copy of the list:
             4ƭ         |     One of these functions, cycling between them:
   ¹                    |       Identity
    ÆlÆeƭ               |       Alternate between log and exp
         ²½ƭ            |       Alternate between square and square root
            İ           |       Reciprocal
               Æm       |    Then take the mean
                       Ḣ| Finally take the first item
Nick Kennedy
la source
Je suis assez mauvais en gelée, mais est-ce que quelque chose de similaire pourrait P*İLfonctionner pour la moyenne géométrique?
mon pronom est monicareinstate
@ Quelqu'un aurait-il besoin de P*Lİ$ne pas économiser des octets? Cela voudrait dire que je pourrais ramener Æmune ligne sans chiffrer d'octets, mais j'aime bien le fait que chacune d'entre elles a actuellement une moyenne arithmétique.
Nick Kennedy
2

Python 3 , 152 octets

from math import*
s=sum
def f(*a):l=len(a);return 2>len({*a})and{*a}or f(s(a)/l,l/s(map(pow,a,l*[-1])),exp(s(map(log,a))/l),(s(map(pow,a,l*[2]))/l)**.5)

Essayez-le en ligne!

Fonction récursive f, convergera vers la précision en virgule flottante. Fonctionne en principe pour toutes les listes de nombres positifs de toute taille, mais est limitée par la récursivité de Python, limitant l’ erreur d’arrondi pour certains cas de test.

Vous pouvez également choisir une précision de 9 décimales:

Python 3 , 169 octets

from math import*
s=sum
def f(*a):l=len(a);return(2>len({round(i,9)for i in a}))*a[0]or f(s(a)/l,l/s(map(pow,a,l*[-1])),exp(s(map(log,a))/l),(s(map(pow,a,l*[2]))/l)**.5)

Essayez-le en ligne!

Jitse
la source
1

C # , 173 octets

double m(int n,params double[]a)=>(n<1?a[0]:m(n-1,a.Sum()/a.Length,Math.Pow(a.Aggregate((t,x)=>t*x),1.0/a.Length),a.Length/a.Sum(x=>1/x),Math.Sqrt(a.Sum(x=>x*x)/a.Length)));

Essayez-le en ligne!

aviateur
la source
2
Cela semble vraiment sur une variable qui doit être passée. En outre, vous devez inclure using Systemet using System.Linqdans votre nombre d'octets, car ils sont nécessaires à l'exécution du programme. Vous pouvez remplacer votre compilateur par le compilateur C # Visual Interactive, qui n'a pas besoin de ces importations. Aussi, 1.0->1d
Incarnation de l'Ignorance
1

Propre , 124 octets

import StdEnv
f=avg o limit o iterate\l=let n=toReal(length l)in[avg l,prod l^(1.0/n),n/sum[1.0/x\\x<-l],avg[x*x\\x<-l]^0.5]

Essayez-le en ligne!

Effectue l'opération jusqu'à ce que le résultat cesse de changer.

Vive la virgule flottante de précision limitée!

Οurous
la source
1

Pyth, 32 octets

h.Wt{H[.OZ@*[email protected]^R2Z2

Essayez-le en ligne ici , ou vérifiez tous les cas de test (sauf deux, voir la note ci-dessous) à la fois ici . Accepte les entrées sous forme de liste.

Il semble y avoir quelques problèmes d'arrondi, car certaines entrées ne convergent pas correctement lorsqu'elles le devraient autrement. En particulier, le scénario de test 0.01 100reste bloqué aux valeurs [6.748305820749738, 6.748305820749738, 6.748305820749739, 6.748305820749738]et le scénario de test 1.61 2.41reste bloqué à[1.9896543776640825, 1.9896543776640825, 1.9896543776640827, 1.9896543776640825] - notez dans les deux cas que la 3ème moyenne (moyenne harmonique) diffère des autres.

Je ne sais pas si ce problème invalide mon inscription, mais je le publie quand même, car cela devrait fonctionner. Si cela n’est pas acceptable, vous pouvez résoudre ce problème en procédant .RRTavant le [, en arrondissant chacun des moyennes à 10 décimales, comme indiqué dans cette suite de tests .

h.Wt{H[.OZ@*[email protected]^R2Z2)Q   Implicit: Q=eval(input())
                                     Trailing )Q inferred
 .W                              Q   Funcitonal while: While condition is true, call inner. Starting value Q
   t{H                               Condition function: current input H
    {H                                 Deduplicate H
   t                                   Discard first value
                                         Empty list is falsey, so while is terminated when means converge
      [.OZ@*[email protected]^R2Z2)    Inner function: current input Z
              JlZ                      Take length of Z, store in J
       .OZ                             (1) Arithmetic mean of Z
           *FZ                         Product of Z
          @   J                        (2) Jth root of the above
                     L Z               Map each element of Z...
                    c 1                ... to its reciprocal
                   s                   Sum the above
                 cJ                    (3) J / the above
                            R Z        Map each element of Z...
                           ^ 2         ... to its square
                         .O            Arithmetic mean of the above
                        @      2       (4) Square root of the above
      [                         )      Wrap results (1), (2), (3), and (4) in a list
                                         This is used as the input for the next iteration of the loop
h                                    Take the first element of the result, implicit print
Sok
la source
Puisque je suis à peu près sûr que le calcul répété ne sautera pas aux valeurs précédentes, vous pouvez le remplacer .Wt{Hpar u-4 octets (et le remplacer Zpar G)
ar4093
1

Japt v2.0a0 -g, 42 38 octets

â ÊÉ?ß[Ux²÷(V=UÊ)¬Ux÷V U×qV V÷Ux!÷1]:U

Il doit y avoir un moyen plus court ... C'est une monstruosité! Sauvegardé 4 octets grâce à @Shaggy!

L'essayer

Incarnation de l'ignorance
la source
38 octets . Mais je suis d'accord, il doit y avoir un chemin plus court!
Shaggy
1

C # (compilateur interactif Visual C #) , 177 octets

double f(double[]g)=>g.All(c=>Math.Abs(c-g[0])<1e-9)?g[0]:f(new[]{g.Sum()/(z=g.Length),Math.Pow(g.Aggregate((a,b)=>a*b),1d/z),z/g.Sum(x=>1/x),Math.Sqrt(g.Sum(x=>x*x)/z)});int z;

Merci à @KevinCruijjsen d'avoir signalé que l'utilisation de la précision en virgule flottante posait problème! Ce serait 163 octets si les doubles étaient parfaitement précis.

Essayez-le en ligne!

Incarnation de l'ignorance
la source
Les deux derniers cas de test donnent une StackOverflowExceptionprécision due à la virgule flottante. Au lieu de c==g[0]vous pourrait faire quelque chose comme Math.Abs(c-g[0])<1e-9. Essayez-le en ligne.
Kevin Cruijssen
@KevinCruijssen Merci, c'est tellement pénible de gérer les nombres à virgule flottante
Incarnation de l'Ignorance
1

Code machine x86 (SIMD 4x float utilisant SSE1 et AVX 128 bits) 94 octets

Code machine x86 (SIMD 4x double utilisant AVX 256 bits) 123 octets

float réussit les tests dans la question, mais avec un seuil de sortie de boucle suffisamment petit pour que cela se produise, il est facile pour lui de rester bloqué dans une boucle infinie avec des entrées aléatoires.

Les instructions SSE1 avec une seule précision empaquetée ont une longueur de 3 octets, mais les instructions SSE2 et AVX simples ont une longueur de 4 octets. (Les instructions Scalar-single, par exemple, sqrtssont également une longueur de 4 octets, c'est pourquoi je l'utilise sqrtpsmême si je ne me soucie que de l'élément low. Ce n'est même pas plus lent que sqrtss sur du matériel moderne). J'ai utilisé AVX pour une destination non destructive pour économiser 2 octets contre movaps + op.
Dans la version double, nous pouvons toujours faire un couple movlhpspour copier des morceaux de 64 bits (car souvent, nous ne nous soucions que de l'élément bas d'une somme horizontale). La somme horizontale d'un vecteur SIMD 256 bits nécessite également un extra vextractf128pour obtenir la moitié haute, par rapport à la stratégie 2x lente mais lente haddpspour float . ledouble version 2x nécessite également des constantes 2x à 8 octets, au lieu de 2x à 4 octets. Globalement, il arrive à près de 4/3 de la taille de lafloat version.

mean(a,b) = mean(a,a,b,b)pour tous les 4 de ces moyens , nous pouvons donc simplement dupliquer l'entrée jusqu'à 4 éléments et ne jamais avoir à mettre en œuvre length = 2. Ainsi, nous pouvons coder en dur une moyenne géométrique telle que 4th-root = sqrt (sqrt), par exemple. Et nous avons seulement besoin d' une constante FP, 4.0.

Nous avons un seul vecteur SIMD sur 4 [a_i, b_i, c_i, d_i]. À partir de cela, nous calculons les 4 moyennes sous forme de scalaires dans des registres séparés, puis nous les remettons ensemble pour la prochaine itération. (Les opérations horizontales sur les vecteurs SIMD ne sont pas pratiques, mais nous devons faire la même chose pour les 4 éléments dans un nombre suffisant de cas. Cela a été équilibré. J'ai commencé avec une version x87 de ce logiciel, mais cela devenait très long et amusant.)

La condition boucle sortie }while(quadratic - harmonic > 4e-5)(ou une plus petite constante double) est basée sur @ de RobinRyder réponse R et la réponse Java Kevin Cruijssen : moyenne quadratique est toujours la plus grande amplitude et moyenne harmonique est toujours la plus petite (abstraction faite des erreurs d' arrondi). Donc, nous pouvons vérifier le delta entre ces deux pour détecter la convergence. Nous renvoyons la moyenne arithmétique en tant que résultat scalaire. C'est généralement entre ces deux et est probablement le moins susceptible d'erreurs d'arrondi.

Version float : appelable comme float meanmean_float_avx(__m128);avec les arguments arg et return en xmm0. (Donc, x86-64 System V, ou Windows x64 vectorcall, mais pas x64 fastcall.) Ou déclarez le type de retour __m128afin que vous puissiez obtenir la moyenne quadratique et harmonique pour le test.

Laisser cela prendre 2 floatarguments distincts dans xmm0 et xmm1 coûterait un octet supplémentaire: il nous faudrait un shufpsavec un imm8 (au lieu de simplement unpcklps xmm0,xmm0) pour mélanger et dupliquer deux entrées.

    40  address                    align 32
    41          code bytes         global meanmean_float_avx
    42                             meanmean_float_avx:
    43 00000000 B9[52000000]           mov      ecx, .arith_mean      ; allows 2-byte call reg, and a base for loading constants
    44 00000005 C4E2791861FC           vbroadcastss  xmm4, [rcx-4]    ; float 4.0
    45                             
    46                                 ;; mean(a,b) = mean(a,b,a,b) for all 4 types of mean
    47                                 ;; so we only ever have to do the length=4 case
    48 0000000B 0F14C0                 unpcklps xmm0,xmm0          ; [b,a] => [b,b,a,a]
    49                             
    50                                 ; do{ ... } while(quadratic - harmonic > threshold);
    51                             .loop:
    52                             ;;; XMM3 = geometric mean: not based on addition.  (Transform to log would be hard.  AVX512ER has exp with 23-bit accuracy, but not log.  vgetexp = floor(lofg2(x)), so that's no good.)
    53                                 ;; sqrt once *first*, making magnitudes closer to 1.0 to reduce rounding error.  Numbers are all positive so this is safe.
    54                                 ;; both sqrts first was better behaved, I think.
    55 0000000E 0F51D8                 sqrtps   xmm3, xmm0                 ; xmm3 = 4th root(x)
    56 00000011 F30F16EB               movshdup xmm5, xmm3                 ; bring odd elements down to even
    57 00000015 0F59EB                 mulps    xmm5, xmm3
    58 00000018 0F12DD                 movhlps  xmm3, xmm5                 ; high half -> low
    59 0000001B 0F59DD                 mulps    xmm3, xmm5                 ; xmm3[0] = hproduct(sqrt(xmm))
    60                             ;    sqrtps   xmm3, xmm3                 ; sqrt(hprod(sqrt)) = 4th root(hprod)
    61                                 ; common final step done after interleaving with quadratic mean
    62                             
    63                             ;;; XMM2 = quadratic mean = max of the means
    64 0000001E C5F859E8               vmulps   xmm5, xmm0,xmm0
    65 00000022 FFD1                   call     rcx                ; arith mean of squares
    66 00000024 0F14EB                 unpcklps xmm5, xmm3         ; [quad^2, geo^2, ?, ?]
    67 00000027 0F51D5                 sqrtps   xmm2, xmm5         ; [quad,   geo,   ?, ?]
    68                             
    69                             ;;; XMM1 = harmonic mean = min of the means
    70 0000002A C5D85EE8               vdivps   xmm5, xmm4, xmm0    ; 4/x
    71 0000002E FFD1                   call     rcx                ; arithmetic mean (under inversion)
    72 00000030 C5D85ECD               vdivps   xmm1, xmm4, xmm5    ; 4/.  (the factor of 4 cancels out)
    73                             
    74                             ;;; XMM5 = arithmetic mean
    75 00000034 0F28E8                 movaps   xmm5, xmm0
    76 00000037 FFD1                   call     rcx
    77                             
    78 00000039 0F14E9                 unpcklps  xmm5, xmm1           ;     [arith, harm, ?,?]
    79 0000003C C5D014C2               vunpcklps xmm0, xmm5,xmm2      ; x = [arith, harm, quad, geo]
    80                             
    81 00000040 0F5CD1                 subps    xmm2, xmm1        ; largest - smallest mean: guaranteed non-negative
    82 00000043 0F2E51F8               ucomiss  xmm2, [rcx-8]     ; quad-harm > convergence_threshold
    83 00000047 73C5                   jae     .loop
    84                             
    85                                 ; return with the arithmetic mean in the low element of xmm0 = scalar return value
    86 00000049 C3                     ret
    87                             
    88                             ;;; "constant pool" between the main function and the helper, like ARM literal pools
    89 0000004A ACC52738           .fpconst_threshold:   dd 4e-5    ; 4.3e-5 is the highest we can go and still pass the main test cases
    90 0000004E 00008040           .fpconst_4:    dd 4.0
    91                             .arith_mean:               ; returns XMM5 = hsum(xmm5)/4.
    92 00000052 C5D37CED               vhaddps   xmm5, xmm5         ; slow but small
    93 00000056 C5D37CED               vhaddps   xmm5, xmm5
    94 0000005A 0F5EEC                 divps     xmm5, xmm4        ; divide before/after summing doesn't matter mathematically or numerically; divisor is a power of 2
    95 0000005D C3                     ret

    96 0000005E 5E000000           .size:      dd $ - meanmean_float_avx
       0x5e = 94 bytes

(Liste NASM créée avec nasm -felf64 mean-mean.asm -l/dev/stdout | cut -b -34,$((34+6))-. Supprimez la partie liste et récupérez le source cut -b 34- > mean-mean.asm)

SIMD La somme horizontale et la division par 4 (c'est-à-dire une moyenne arithmétique) sont implémentées dans une fonction distincte que nous call(avec un pointeur de fonction pour amortir le coût de l'adresse). Avec 4/xavant / après, ou x^2avant et après après, nous obtenons la moyenne harmonique et la moyenne quadratique. (Il était pénible d’écrire ces divinstructions au lieu de les multiplier par un exactement représentable 0.25.)

La moyenne géométrique est mise en œuvre séparément avec multiplie et sqrt chaîné. Ou avec un premier carré pour réduire la magnitude des exposants et peut-être aider la précision numérique. le journal n'est pas disponible, uniquement floor(log2(x))via AVX512vgetexpps/pd . Exp est en quelque sorte disponible via AVX512ER (Xeon Phi uniquement), mais avec une précision de 2 ^ -23 seulement.

Le mélange d'instructions AVX 128 bits et de SSE hérité n'est pas un problème de performances. Le mélange d’AVX 256 bits et de SSE hérité peut s’effectuer sur Haswell, mais sur Skylake, cela crée potentiellement une fausse dépendance potentielle pour les instructions SSE. Je pense que ma doubleversion évite les chaînes dep inutiles transportées en boucle, ainsi que les goulots d'étranglement sur la latence / le débit de div / sqrt.

Version double:

   108                             global meanmean_double_avx
   109                             meanmean_double_avx:
   110 00000080 B9[E8000000]           mov      ecx, .arith_mean
   111 00000085 C4E27D1961F8           vbroadcastsd  ymm4, [rcx-8]    ; float 4.0
   112                             
   113                                 ;; mean(a,b) = mean(a,b,a,b) for all 4 types of mean
   114                                 ;; so we only ever have to do the length=4 case
   115 0000008B C4E37D18C001           vinsertf128   ymm0, ymm0, xmm0, 1       ; [b,a] => [b,a,b,a]
   116                             
   117                             .loop:
   118                             ;;; XMM3 = geometric mean: not based on addition.
   119 00000091 C5FD51D8               vsqrtpd      ymm3, ymm0     ; sqrt first to get magnitude closer to 1.0 for better(?) numerical precision
   120 00000095 C4E37D19DD01           vextractf128 xmm5, ymm3, 1           ; extract high lane
   121 0000009B C5D159EB               vmulpd       xmm5, xmm3
   122 0000009F 0F12DD                 movhlps      xmm3, xmm5              ; extract high half
   123 000000A2 F20F59DD               mulsd        xmm3, xmm5              ; xmm3 = hproduct(sqrt(xmm0))
   124                                ; sqrtsd       xmm3, xmm3             ; xmm3 = 4th root = geomean(xmm0)   ;deferred until quadratic
   125                             
   126                             ;;; XMM2 = quadratic mean = max of the means
   127 000000A6 C5FD59E8               vmulpd   ymm5, ymm0,ymm0
   128 000000AA FFD1                   call     rcx                ; arith mean of squares
   129 000000AC 0F16EB                 movlhps  xmm5, xmm3         ; [quad^2, geo^2]
   130 000000AF 660F51D5               sqrtpd   xmm2, xmm5         ; [quad  , geo]
   131                             
   132                             ;;; XMM1 = harmonic mean = min of the means
   133 000000B3 C5DD5EE8               vdivpd   ymm5, ymm4, ymm0    ; 4/x
   134 000000B7 FFD1                   call     rcx                 ; arithmetic mean under inversion
   135 000000B9 C5DB5ECD               vdivsd   xmm1, xmm4, xmm5    ; 4/.  (the factor of 4 cancels out)
   136                             
   137                             ;;; XMM5 = arithmetic mean
   138 000000BD C5FC28E8               vmovaps  ymm5, ymm0
   139 000000C1 FFD1                   call     rcx
   140                             
   141 000000C3 0F16E9                 movlhps     xmm5, xmm1            ;     [arith, harm]
   142 000000C6 C4E35518C201           vinsertf128 ymm0, ymm5, xmm2, 1   ; x = [arith, harm, quad, geo]
   143                             
   144 000000CC C5EB5CD1               vsubsd   xmm2, xmm1               ; largest - smallest mean: guaranteed non-negative
   145 000000D0 660F2E51F0             ucomisd  xmm2, [rcx-16]           ; quad - harm > threshold
   146 000000D5 77BA                   ja      .loop
   147                             
   148                                 ; vzeroupper ; not needed for correctness, only performance
   149                                 ; return with the arithmetic mean in the low element of xmm0 = scalar return value
   150 000000D7 C3                     ret
   151                             
   152                             ; "literal pool" between the function
   153 000000D8 95D626E80B2E113E   .fpconst_threshold:   dq 1e-9
   154 000000E0 0000000000001040   .fpconst_4:    dq 4.0            ; TODO: golf these zeros?  vpbroadcastb and convert?
   155                             .arith_mean:                     ; returns YMM5 = hsum(ymm5)/4.
   156 000000E8 C4E37D19EF01           vextractf128 xmm7, ymm5, 1
   157 000000EE C5D158EF               vaddpd       xmm5, xmm7
   158 000000F2 C5D17CED               vhaddpd      xmm5, xmm5      ; slow but small
   159 000000F6 C5D35EEC               vdivsd     xmm5, xmm4        ; only low element matters
   160 000000FA C3                     ret

   161 000000FB 7B000000           .size:      dd $ - meanmean_double_avx

    0x7b = 123 bytes

Harnais de test C

#include <immintrin.h>
#include <stdio.h>
#include <math.h>

static const struct ab_avg {
    double a,b;
    double mean;
} testcases[] = {
    {1, 1, 1},
    {1, 2, 1.45568889},
    {100, 200, 145.568889},
    {2.71, 3.14, 2.92103713},
    {0.57, 1.78, 1.0848205},
    {1.61, 2.41, 1.98965438},
    {0.01, 100, 6.7483058},
};

// see asm comments for order of  arith, harm, quad, geo
__m128 meanmean_float_avx(__m128);       // or float ...
__m256d meanmean_double_avx(__m128d);    // or double ...
int main(void) {
    int len = sizeof(testcases) / sizeof(testcases[0]);
    for(int i=0 ; i<len ; i++) {
        const struct ab_avg *p = &testcases[i];
#if 1
        __m128 arg = _mm_set_ps(0,0, p->b, p->a);
        double res = meanmean_float_avx(arg)[0];
#else
        __m128d arg = _mm_loadu_pd(&p->a);
        double res = meanmean_double_avx(arg)[0];
#endif
        double allowed_diff = (p->b - p->a) / 100000.0;
        double delta = fabs(p->mean - res);
        if (delta > 1e-3 || delta > allowed_diff) {
            printf("%f %f => %.9f but we got %.9f.  delta = %g allowed=%g\n",
                   p->a, p->b, p->mean, res, p->mean - res, allowed_diff);
        }
    }



    while(1) {
        double a = drand48(), b = drand48();  // range= [0..1)
        if (a>b) {
            double tmp=a;
            a=b;
            b=tmp; // sorted
        }
//      a *= 0.00000001;
//      b *= 123156;
        // a += 1<<11;  b += (1<<12)+1;  // float version gets stuck inflooping on 2048.04, 4097.18 at fpthreshold = 4e-5

        // a *= 1<<11 ; b *= 1<<11;   // scaling to large magnitude makes sum of squares loses more precision
        //a += 1<<11; b+= 1<<11;   // adding to large magnitude is hard for everything, catastrophic cancellation
#if 1
        printf("testing float %g, %g\n", a, b);
        __m128 arg = _mm_set_ps(0,0, b, a);
        __m128 res = meanmean_float_avx(arg);
        double quad = res[2], harm = res[1];  // same order as double... for now
#else
        printf("testing double %g, %g\n", a, b);
        __m128d arg = _mm_set_pd(b, a);
        __m256d res = meanmean_double_avx(arg);
        double quad = res[2], harm = res[1];
#endif
        double delta = fabs(quad - harm);
        double allowed_diff = (b - a) / 100000.0; // calculated in double even for the float case.
        // TODO: use the double res as a reference for float res
        // instead of just checking quadratic vs. harmonic mean

        if (delta > 1e-3 || delta > allowed_diff) {
            printf("%g %g we got q=%g, h=%g, a=%g.  delta = %g,  allowed=%g\n",
                   a, b, quad, harm, res[0], quad-harm, allowed_diff);
        }
    }

}

Construire avec:

nasm -felf64 mean-mean.asm &&
gcc -no-pie -fno-pie -g -O2 -march=native mean-mean.c mean-mean.o

Évidemment, vous avez besoin d’un processeur prenant en charge AVX ou d’un émulateur comme Intel SDE. Pour compiler sur un hôte sans support natif d’AVX, utilisez -march=sandybridgeou-mavx

Exécuter: réussit les cas de test codés en dur, mais pour la version flottante, les cas de test aléatoires échouent souvent au (b-a)/10000seuil défini dans la question.

$ ./a.out
 (note: empty output before the first "testing float" means clean pass on the constant test cases)
testing float 3.90799e-14, 0.000985395
3.90799e-14 0.000985395 we got q=3.20062e-10, h=3.58723e-05, a=2.50934e-05.  delta = -3.5872e-05,  allowed=9.85395e-09
testing float 0.041631, 0.176643
testing float 0.0913306, 0.364602
testing float 0.0922976, 0.487217
testing float 0.454433, 0.52675
0.454433 0.52675 we got q=0.48992, h=0.489927, a=0.489925.  delta = -6.79493e-06,  allowed=7.23169e-07
testing float 0.233178, 0.831292
testing float 0.56806, 0.931731
testing float 0.0508319, 0.556094
testing float 0.0189148, 0.767051
0.0189148 0.767051 we got q=0.210471, h=0.210484, a=0.21048.  delta = -1.37389e-05,  allowed=7.48136e-06
testing float 0.25236, 0.298197
0.25236 0.298197 we got q=0.274796, h=0.274803, a=0.274801.  delta = -6.19888e-06,  allowed=4.58374e-07
testing float 0.531557, 0.875981
testing float 0.515431, 0.920261
testing float 0.18842, 0.810429
testing float 0.570614, 0.886314
testing float 0.0767746, 0.815274
testing float 0.118352, 0.984891
0.118352 0.984891 we got q=0.427845, h=0.427872, a=0.427863.  delta = -2.66135e-05,  allowed=8.66539e-06
testing float 0.784484, 0.893906
0.784484 0.893906 we got q=0.838297, h=0.838304, a=0.838302.  delta = -7.09295e-06,  allowed=1.09422e-06

Les erreurs de PF suffisent pour que quad-harm soit inférieur à zéro pour certaines entrées.

Ou sans a += 1<<11; b += (1<<12)+1;commentaire:

testing float 2048, 4097
testing float 2048.04, 4097.18
^C  (stuck in an infinite loop).

Aucun de ces problèmes n'arrive avec double. Mettez en commentaire l' printfavant chaque test pour vérifier que la sortie est vide (rien du if(delta too high)bloc).

TODO: utilisez la doubleversion comme référence pour la floatversion, au lieu de simplement regarder comment elles convergent avec quad-harm.

Peter Cordes
la source
1

Javascript - 186 octets

Prend les entrées sous forme de tableau de nombres. Utilise les transformations moyennes dans la réponse de J42161217 pour raccourcir le code.

Essayez-le en ligne

f=(v,l=[m=(w,k=0)=>w.map(x=>k+=x)&&k/w.length,w=>1/m(w.map(x=>1/x)),w=>Math.E**m(w.map(x=>Math.log(x))),w=>m(w.map(x=>x**2))**.5].map(x=>x(v)).sort((a,b)=>a-b))=>l[3]-l[0]>1e-5?f(l):l[0]

Explication

f = (
  v,
  l=[
    m=(w,k=0)=>w.map(x=>k+=x)&&k/w.length,  // m = w => arithmetic mean of values in w
    w=>1/m(w.map(x=>1/x)),                  // w => harmonic mean of values in w   
    w=>Math.E**m(w.map(x=>Math.log(x))),    // w => geometric mean of values in w   
    w=>m(w.map(x=>x**2))**.5                // w => quadratic mean of values in w   
  ].map(x=>x(v))                            // get array of each mean using input v, stored in l
  .sort((a,b)=>a-b)                         // sort the outputs
) =>
  l[3] - l[0] > 1e-5 ?                      // is the difference between the largest
                                            // and smallest means > 1/100000?
    f(l) :                                  // if yes, get the mean mean of the means
    l[0]                                    // if no, arbitrarily return the smallest value
                                            // as close enough
marchand
la source
Je pensais que j'allais être intelligent et mettre en œuvre la relation avec les logarithmes, mais il semblerait que vous et J42161217 y soient arrivés les premiers!
Pureferret le
@Pureferret Je ne prends aucun crédit pour cela, je l'ai volé de manière flagrante: D
asgallant
vous l'avez cependant écrit en JavaScript!
Pureferret le
1
C'était la partie facile. C'était difficile de jouer au golf.
asgallant le
1
Le TIL n'a pas été configuré correctement. J'ai ajouté un lien TIL à la réponse.
asgallant
0

SNOBOL4 (CSNOBOL4) , 296 octets

	X =INPUT
	Y =INPUT
	A =(X + Y) / 2.
	P =X * Y
	G =P ^ 0.5
	H =P / A
	Q =(2 * (A ^ 2) - P) ^ 0.5
O	OUTPUT =EQ(Q,A) Q	:S(END)
	M =A
	N =G
	O =H
	P =Q
	A =(M + N + O + P) / 4
	G =(M * N * O * P) ^ 0.25
	H =4 / (1 / M + 1 / N + 1 / O + 1 / P)
	Q =((M ^ 2 + N ^ 2 + O ^ 2 + P ^ 2) / 4) ^ 0.5	:(O)
END

Essayez-le en ligne!

Mise en œuvre simple. Utilise un truc de ma réponse à une question liée au golf un peu plus.

Giuseppe
la source