Berechnen Sie den Mittelwert zweier Zahlen

41

Haftungsausschluss: Der Mittelwert wird von mir gebildet

Definiere das arithmetische Mittel von Zahlen als Definiere das geometrische Mittel von Zahlen als Definieren Sie das harmonische Mittel von Zahlen als Definiere den quadratischen Mittelwert von Zahlen als Der Mittelwert ( ) ist wie folgt definiert: Definiere vier Folgen ( ) alsn

M1(x1,...,xn)=x1+x2+...+xnn
n
M0(x1,...,xn)=x1x2...xnn
n
M1(x1,...,xn)=n1x2+1x2+...+1xn
n
M2(x1,...,xn)=x12+x22+...+xn2n
MMak,bk,ck,dk
a0=M1(x1,...,xn),b0=M0(x1,...,xn),c0=M1(x1,...,xn),d0=M2(x1,...,xn),ak+1=M1(ak,bk,ck,dk),bk+1=M0(ak,bk,ck,dk),ck+1=M1(ak,bk,ck,dk),dk+1=M2(ak,bk,ck,dk)
M M ( x 1 , x 2 , . . . , x n ) Alle vier Sequenzen konvergieren zu die gleiche Zahl, .MM(x1,x2,...,xn)

Beispiel

Der Mittelwert von 1 und 2 wird wie folgt berechnet: Beginnen Sie mit Dann ist Die weitere Berechnung der Sequenzen sollte klar sein. Es ist ersichtlich, dass sie zu derselben Zahl, ungefähr , .

a0=(1+2)/2=1.5,b0=12=21.4142,c0=211+12=431.3333,d0=12+222=521.5811.
a1=1.5+1.4142+1.3333+1.581141.4571,b1=1.51.41421.33331.581141.4542,c1=411.5+11.4142+11.3333+11.58111.4512,d1=1.52+1.41422+1.33332+1.5811241.4601.
1.45568889

Herausforderung

Berechnen Sie bei zwei positiven reellen Zahlen und ( ) ihren Mittelwert .aba<bMM(a,b)

Testfälle

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

Anmerkungen

  • Ihr Programm ist gültig, wenn die Differenz zwischen seiner Ausgabe und der korrekten Ausgabe nicht größer als 1/100000 des absoluten Wertes der Differenz zwischen den Eingabezahlen ist.
  • Die Ausgabe sollte eine einzelne Zahl sein.

Das ist , also gewinnt der kürzeste Code!

jemand
quelle
Sandbox-Link
jemand
11
Wie genau sollen wir sein?
Verkörperung der Ignoranz
2
Eng verwandt
Giuseppe
1
Können wir annehmen, dass die erste Eingabe immer kleiner als die zweite ist, wie in all Ihren Testfällen? (Wenn nicht, werde ich meine Java-Antwort
zurücksetzen

Antworten:

14

Wolfram Language (Mathematica) , 52 Byte

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

Probieren Sie es online!

In meinem ersten Ansatz habe ich diese eingebauten
Mean GeometricMean HarmonicMeanund verwendetRootMeanSquare

Hier sind einige Ersetzungen zum Speichern von Bytes

HarmonicMean-> 1/Mean[1/x] von @Robin Ryder (3 Bytes gespeichert)
GeometricMean-> E^Mean@Log@xvon @A. Rex (2 Bytes gespeichert)
RootMeanSquare-> Mean[x^2]^.5von @A. Rex (4 Bytes gespeichert)

schließlich können wir zuweisen Meanzu M(wie von @ovs vorgeschlagen) und speichern 5 mehr Bytes

J42161217
quelle
Sparen Sie 2 Bytes, indem Sie GeometricMean neu codieren
Robin Ryder
@ RobinRyder Ich glaube du meinst Harmonic .. nett!
J42161217
1
Speichern Sie 8 weitere Bytes :#//.x_:>N@{Mean@x,E^Mean@Log@x,1/Mean[1/x],Mean[x^2]^.5}&
A. Rex
@ovs bearbeitet .....
J42161217
10

R 70 69 67 Bytes

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

Probieren Sie es online!

-1 Byte bei besserer Konditionierung.
-2 Bytes durch Umschalten auf Basis 2.

Wie bei einigen anderen Antworten wird der Ausdruck des geometrischen Mittels als arithmetisches Mittel auf der logarithmischen Skala (hier in Basis 2) verwendet:

M0(x1,,xn)=2M1(log2x1,,log2xn).

Es wird auch die Tatsache verwendet , dass , dh . Die Bedingung ist daher äquivalent zu , was ich in der while-Schleife verwende; Dies wird erreicht, indem die Syntax missbraucht wird, die nur das erste Element berücksichtigt, wenn die Bedingung ein Vektor ist, daher die Reihenfolge, in der die Mittel gespeichert sind. (Beachten Sie, dass wir stattdessen auch verwenden da dies das Minimum der vier ist, aber wir könnten in der Bedingung nicht oder .)k,dkakbkckdk=max(ak,bk,ck,dk)ak=bk=ck=dkdk=M1(ak,bk,ck,dk)whileckakbk

Wenn wir die while-Schleife verlassen, xist dies ein konstanter Vektor. Das Finale ?xberechnet seinen Mittelwert, um ihn auf einen Skalar zu reduzieren.

Robin Ryder
quelle
1
Sollte es nicht anstelle von ? l o g x nlnxnlogxn
Tau
@Tau Ja, ich habe den natürlichen Logarithmus durch , was in R der Standard ist. Wie auch immer, ich habe ihn jetzt auf Basis 2 Logarithmus für -2 Bytes geändert. log
Robin Ryder
6

J , 34 Bytes

(31 als Ausdruck ohne Zuweisung zur Variablen f)

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

Probieren Sie es online!

Für Funktionen , aund b, a &.: b( "a unter b" ( verwandte Herausforderung )) entspricht (b inv) a b- gelten b, dann a, dann Inverse von b. In diesem Fall ist der geometrische / harmonische / quadratische Mittelwert der arithmetische Mittelwert "unter" Logarithmus, Inversion bzw. Quadrat.

user202729
quelle
5

TI-BASIC, 42 35 34 Bytes

-1 Byte dank @SolomonUcko

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

Die Eingabe ist eine Liste von zwei Ganzzahlen in Ans.
Die Ausgabe wird in gespeichert Ansund nach Abschluss des Programms automatisch ausgedruckt.

Formeln, die für geometrische, harmonische und quadratische Mittel verwendet werden, basieren auf der Erklärung von user202729 .

Beispiel:

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

Erläuterung:
(Zeilenumbrüche wurden zur Verdeutlichung hinzugefügt. Sie erscheinen NICHT im 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

Anmerkungen:

TI-BASIC ist eine Token-Sprache. Die Anzahl der Zeichen entspricht nicht der Anzahl der Bytes.

e^(ist dieses Ein-Byte-Token.

^-1wird für dieses Ein-Byte-Token verwendet.
Ich habe mich ^-1stattdessen für das Schreiben entschieden, weil das Token so aussieht, ֿ¹als ob es sich in einem Codeblock befindet.

√(ist dieses Ein-Byte-Token.

ΔList(ist dieses Zwei-Byte-Token.

Tau
quelle
Ich denke, Sie können eine Klammer retten, indem Sie den geometrischen Mittelwert zuletzt setzen.
Solomon Ucko
@SolomonUcko ah, danke fürs mitbekommen! Ich habe das vorher nicht bedacht.
Tau
max(DeltaList(Ans-> variance(Ans.
Lirtosiast
5

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

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 Bytes dank @PeterCordes .
-15 weitere Bytes dank @PeterCordes , inspiriert von der R-Antwort von @RobinRyder .
+4 Byte, da ich davon ausgegangen bin, dass die Eingänge vorbestellt sind.
-27 Bytes dank @ OlivierGrégoire .

Probieren Sie es online aus.

Erläuterung:

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
quelle
In C könnten Sie f+=Math.abs(d-D)<1e-9;eine implizite Konvertierung von einem booleschen Vergleichsergebnis in eine 0/1-Ganzzahl erhalten und dann double. Hat Java eine kompakte Syntax dafür? Oder ist es möglich, zu f+=Math.abs(d-D)überprüfen, ob die Summe der absoluten Differenzen klein genug ist?
Peter Cordes
1
Yup arbeitet für Ihre Testfälle f>1e-8als Schleifenbedingung: 229 Bytes. 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];}. Mit 1e-9läuft es langsamer (ungefähr doppelt so lange wie die CPU-Zeit) und muss mehr Iterationen ausführen, um im Wesentlichen 4 * d-Dso klein zu machen. Mit 1e-7ist es ungefähr die gleiche Geschwindigkeit wie 1e-8. Mit 1e-6unterscheiden sich einige der nachfolgenden Ziffern für einen Fall.
Peter Cordes
1
@ RobinRyders Antwort weist darauf hin, dass der quadratische Mittelwert immer der größte und der harmonische immer der kleinste ist. Vielleicht können Sie also fganz aufgeben und nur nachprüfen a[3]-a[2]<4e-9.
Peter Cordes
1
@PeterCordes l==2||du meinst (Golfen l<3|). Aber ja, guter Punkt; Ich habe es hinzugefügt. :)
Kevin Cruijssen
2
180 Bytes durch Aggregation aggregierbarer Reduktionen.
Olivier Grégoire
3

Kohle , 40 Bytes

W‹⌊θ⌈θ≔⟦∕ΣθLθXΠθ∕¹Lθ∕LθΣ∕¹θ₂∕ΣXθ²Lθ⟧θI⊟θ

Probieren Sie es online! Link ist eine ausführliche Version des Codes. Nimmt die Eingabe als Array von Zahlen. Erläuterung:

W‹⌊θ⌈θ

Wiederholen, während das Array verschiedene Werte enthält ...

≔⟦....⟧θ

... ersetzen Sie das Array durch eine Liste von Werten:

∕ΣθLθ

... die meine ...

XΠθ∕¹Lθ

... das geometrische Mittel ...

∕LθΣ∕¹θ

... das harmonische Mittel ...

₂∕ΣXθ²Lθ

... und das quadratische Mittel.

I⊟θ

Konvertieren Sie ein Element des Arrays in einen String und drucken Sie es implizit aus.

Neil
quelle
3

PowerShell , 182 180 183 Byte

$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]}

Probieren Sie es online!

Andrei Odegov
quelle
171 Bytes
mazzy
@mazzy, beeindruckend :)
Andrei Odegov
3

05AB1E , 26 24 23 Bytes

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

Probieren Sie es online aus oder sehen Sie sich die Schritte aller Testfälle an .

-1 Byte dank @Grimy .

23-Byte-Alternative für geometrisches Mittel:

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

Probieren Sie es online aus oder sehen Sie sich die Schritte aller Testfälle an .

Erläuterung:

Δ         # 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
quelle
23:Δ©P®gzm®ÅA®zÅAz®nÅAt)}н
Grimmy
@ Grimy Danke! Ich kann nicht glauben, dass ich nicht daran gedacht hatte, die Länge Yfür 2/4 zu verwenden. :)
Kevin Cruijssen
1
Weitere 23 , dass betters zeigt die Ähnlichkeit der geometrische Mittel zu den anderen: Δ©ÅA®.²ÅAo®zÅAz®nÅAt)}н. Leider sieht es nicht so aus, als könnten wir all diese Dinge umgestalten ÅA.
Grimmy
@ Grimy Oh, ich mag diese zweite Version. :) EDIT: Ups .. danke, dass du meinen Fehler in der Erklärung bemerkt
hast
Ich programmiere in 05ab1e nicht sehr gut, aber können Sie die Summen berechnen und später alle durch die Länge dividieren?
Jemand
2

Jelly , 25 24 Bytes

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

Probieren Sie es online!

Erläuterung

                    µÐ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
quelle
Ich bin ziemlich schlecht in Jelly, aber könnte etwas Ähnliches P*İLfür den geometrischen Mittelwert funktionieren?
jemand
@jemand müsste es sein P*Lİ$, würde also keine Bytes sparen. Es würde bedeuten, dass ich Æmeine Zeile zurücksetzen könnte, ohne Byte zu kosten, aber ich mag die Tatsache, dass jeder momentan ein arithmetisches Mittel in seinem Kern hat.
Nick Kennedy
2

Python 3 , 152 Bytes

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)

Probieren Sie es online!

Rekursive Funktion f, konvergiert zur Gleitkommapräzision. Funktioniert im Prinzip für alle Listen mit positiven Zahlen beliebiger Größe, ist jedoch durch die Python-Rekursion begrenzt und begrenzt einen Rundungsfehler für einige Testfälle.


Alternativ können Sie sich für eine Genauigkeit von 9 Dezimalstellen entscheiden:

Python 3 , 169 Bytes

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)

Probieren Sie es online!

Jitse
quelle
1

C # 173 Bytes

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)));

Probieren Sie es online!

Flieger
quelle
2
Dies scheint wirklich eine Variable zu sein, die übergeben werden muss. Außerdem müssen Sie using Systemund using System.Linqin Ihre Byteanzahl aufnehmen, da diese für die Ausführung des Programms erforderlich sind. Sie können Ihren Compiler in den C # Visual Interactive Compiler ändern, für den diese Importe nicht erforderlich sind. Auch 1.0->1d
Verkörperung der Ignoranz
1

Sauber , 124 Bytes

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]

Probieren Sie es online!

Führt den Vorgang aus, bis sich das Ergebnis nicht mehr ändert.

Hurra für Gleitkommazahlen mit begrenzter Präzision!

Οurous
quelle
1

Pyth, 32 Bytes

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

Versuchen Sie es online hier oder überprüfen alle die Testfälle (Balken zwei, siehe Anmerkung unten) auf einmal hier . Akzeptiert Eingaben als Liste.

Es scheint einige Rundungsprobleme zu geben, da bestimmte Eingaben nicht richtig konvergieren, wenn dies ansonsten der Fall sein sollte. Insbesondere 0.01 100bleibt der Testfall bei Werten stecken [6.748305820749738, 6.748305820749738, 6.748305820749739, 6.748305820749738], und der Testfall 1.61 2.41bleibt stecken [1.9896543776640825, 1.9896543776640825, 1.9896543776640827, 1.9896543776640825]- beachten Sie in beiden Fällen, dass sich der 3. Mittelwert (harmonischer Mittelwert) von den anderen unterscheidet.

Ich bin nicht sicher, ob dieses Problem meine Eingabe ungültig macht, aber ich poste es trotzdem, wie es funktionieren sollte . Wenn dies nicht akzeptabel ist, kann es behoben werden, indem .RRTvor dem angegeben [wird, dass jedes Mittel auf 10 Dezimalstellen gerundet werden soll, wie in dieser Testsuite dargestellt .

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
quelle
Da ich ziemlich sicher , dass die wiederholte Berechnung nicht zu früheren Werten herumspringen, könnten Sie ersetzen .Wt{Hmit ufür -4 Bytes (und Änderung Zan G)
ar4093
1

Japt v2.0a0 -g, 42 38 Bytes

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

Es muss einen kürzeren Weg geben ... Das ist eine Monstrosität! 4 Bytes gespart dank @Shaggy!

Versuch es

Verkörperung der Ignoranz
quelle
38 Bytes . Aber ich stimme zu, es muss einen kürzeren Weg geben!
Shaggy
1

C # (Visual C # Interactive Compiler) , 177 Byte

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;

Vielen Dank an @KevinCruijjsen für den Hinweis, dass die Verwendung der Gleitkommapräzision Probleme verursacht hat! Wäre 163 Bytes, wenn Doppelte absolut präzise wären.

Probieren Sie es online!

Verkörperung der Ignoranz
quelle
Die letzten beiden Testfälle ergeben StackOverflowExceptionaufgrund der Gleitkommapräzision eine. Anstelle von c==g[0]dir könntest du sowas machen Math.Abs(c-g[0])<1e-9. Probieren Sie es online aus.
Kevin Cruijssen
@ KevinCruijssen Vielen Dank, es ist so ein Schmerz im Umgang mit Gleitkommazahlen
Verkörperung der Ignoranz
1

x86-Maschinencode (SIMD 4x Float mit 128-Bit-SSE1 und AVX) 94 Byte

x86-Maschinencode (SIMD 4x Double mit 256-Bit-AVX) 123 Byte

floatBesteht die Testfälle in der Frage, aber mit einer Loop-Exit-Schwelle, die klein genug ist, um dies zu ermöglichen, kann es leicht passieren, dass der Benutzer in einer Endlosschleife mit zufälligen Eingaben steckt.

SSE1-Befehle mit einfacher Genauigkeit sind 3 Byte lang, SSE2- und einfache AVX-Befehle 4 Byte. (Skalar-Einzelanweisungen wie sqrtsssind auch 4 Byte lang, weshalb ich sie verwende sqrtps, obwohl mir nur das niedrige Element am Herzen liegt. Es ist nicht einmal langsamer als sqrtss auf moderner Hardware). Ich habe AVX als zerstörungsfreies Ziel verwendet, um 2 Bytes gegenüber movaps + op zu sparen.
In der doppelten Version können wir immer noch ein paar movlhps64-Bit-Blöcke kopieren (da uns oft nur das niedrige Element einer horizontalen Summe wichtig ist). Die horizontale Summe eines 256-Bit-SIMD-Vektors erfordert auch ein Extra vextractf128, um die hohe Hälfte zu erhalten, im Vergleich zu der langsamen, aber kleinen 2x- haddpsStrategie für das Floaten . DasdoubleVersion benötigt auch 2x 8-Byte-Konstanten anstelle von 2x 4-Byte. Insgesamt kommt es bei knapp 4/3 der Größe der floatVersion heraus.

mean(a,b) = mean(a,a,b,b)für alle 4 dieser Mittel können wir also einfach die Eingabe auf 4 Elemente duplizieren und müssen nie length = 2 implementieren. So können wir zum Beispiel den geometrischen Mittelwert als 4th-root = sqrt (sqrt) fest codieren. Und wir brauchen nur eine FP - Konstante 4.0.

Wir haben einen einzigen SIMD-Vektor von allen 4 [a_i, b_i, c_i, d_i]. Daraus berechnen wir die 4 Mittelwerte als Skalare in separaten Registern und mischen sie für die nächste Iteration wieder zusammen. (Horizontale Operationen an SIMD-Vektoren sind unpraktisch, aber wir müssen das Gleiche für alle 4 Elemente in genügend Fällen tun, damit es ausgeglichen wird. Ich habe mit einer x87-Version davon begonnen, aber es wurde sehr lang und machte keinen Spaß.)

Die Schleifenausgangsbedingung }while(quadratic - harmonic > 4e-5)(oder eine kleinere Konstante für double) basiert auf @ RobinRyders R-Antwort und Kevin Cruijssens Java-Antwort : Der quadratische Mittelwert ist immer die größte und der harmonische Mittelwert ist immer die kleinste (ohne Berücksichtigung von Rundungsfehlern). So können wir das Delta zwischen diesen beiden prüfen, um Konvergenz zu erkennen. Wir geben das arithmetische Mittel als skalares Ergebnis zurück. Sie liegt normalerweise zwischen diesen beiden Werten und ist wahrscheinlich am wenigsten anfällig für Rundungsfehler.

Float-Version : Aufrufbar wie float meanmean_float_avx(__m128);bei arg und Rückgabewert in xmm0. (Also x86-64 System V oder Windows x64 vectorcall, aber nicht x64 fastcall.) Oder deklarieren Sie den Rückgabetyp so, __m128dass Sie den quadratischen und harmonischen Mittelwert zum Testen erhalten.

floatEs würde 1 zusätzliches Byte kosten, wenn dies 2 separate Argumente in xmm0 und xmm1 annehmen würde: Wir müssten ein shufpsmit einem imm8 (anstatt nur unpcklps xmm0,xmm0) zusammenmischen und 2 Eingaben duplizieren.

    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

(NASM-Auflistung erstellt mit nasm -felf64 mean-mean.asm -l/dev/stdout | cut -b -34,$((34+6))-. Entfernen Sie den Auflistungsteil und stellen Sie die Quelle mit wieder her. cut -b 34- > mean-mean.asm)

SIMD horizontale Summe und Division durch 4 (dh arithmetisches Mittel) wird in einer separaten Funktion implementiert, die wir call(mit einem Funktionszeiger zum Amortisieren der Kosten der Adresse) verwenden. Mit 4/xvorher / nachher oder x^2vorher und nachher erhalten wir den harmonischen Mittelwert und den quadratischen Mittelwert. (Es war schmerzhaft, diese divAnweisungen zu schreiben , anstatt sie mit einer genau darstellbaren Zahl zu multiplizieren 0.25.)

Der geometrische Mittelwert wird separat mit multiplizieren und verkettetem sqrt implementiert. Oder mit einem Quadrat zuerst, um die Größe des Exponenten zu verringern und möglicherweise die numerische Genauigkeit zu verbessern. log ist nicht verfügbar, nur floor(log2(x))über AVX512 vgetexpps/pd. Exp ist irgendwie über AVX512ER verfügbar (nur Xeon Phi), aber nur mit einer Genauigkeit von 2 ^ -23.

Das Mischen von 128-Bit-AVX-Anweisungen und Legacy-SSE ist kein Leistungsproblem. Das Mischen von 256-Bit-AVX mit Legacy-SSE kann in Haswell erfolgen, aber in Skylake kann dies möglicherweise zu einer falschen Abhängigkeit für SSE-Anweisungen führen. Ich denke, meine doubleVersion vermeidet unnötige Loop-Carraged-Dep-Ketten und Engpässe bei der div / sqrt-Latenz / dem Durchsatz.

Doppelversion:

   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

C-Testkabel

#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);
        }
    }

}

Bauen mit:

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

Natürlich benötigen Sie eine CPU mit AVX-Unterstützung oder einen Emulator wie Intel SDE. Verwenden Sie -march=sandybridgeoder, um auf einem Host ohne native AVX-Unterstützung zu kompilieren-mavx

Ausführen: Besteht die hartcodierten Testfälle, aber bei der Float-Version (b-a)/10000überschreiten zufällige Testfälle häufig den in der Frage festgelegten Schwellenwert.

$ ./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

FP-Fehler reichen aus, damit der Quad-Harm bei einigen Eingängen kleiner als Null ist.

Oder mit a += 1<<11; b += (1<<12)+1;unkommentierten:

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

Keines dieser Probleme passiert mit double. Kommentieren Sie printfvor jedem Test aus, um festzustellen, dass der Ausgang leer ist (nichts aus dem if(delta too high)Block).

TODO: Verwenden Sie die doubleVersion als Referenz für die floatVersion, anstatt nur zu betrachten, wie sie mit Quad-Harm konvergieren.

Peter Cordes
quelle
1

Javascript - 186 Bytes

Nimmt die Eingabe als Array von Zahlen. Verwendet die mittleren Transformationen in der Antwort von J42161217 , um den Code zu verkürzen.

Probieren Sie es online

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]

Erläuterung

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
asgallant
quelle
Ich dachte, ich würde klug sein und die Beziehung mit Logarithmen implementieren, aber es sieht so aus, als wären Sie und J42161217 zuerst da!
Pureferret
@ Pureferret Ich nehme keine Anerkennung dafür, ich habe es
krass
Du hast es aber in JavaScript geschrieben!
Pureferret
1
Das war der einfache Teil. Golfen war schwer.
Asgallant
1
Die TIL wurde nicht richtig konfiguriert. Ich habe der Antwort einen TIL-Link hinzugefügt.
Asgallant
0

SNOBOL4 (CSNOBOL4) , 296 Bytes

	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

Probieren Sie es online!

Einfache Implementierung. Benutzt einen Trick von meiner Antwort auf eine verwandte Frage, um ein bisschen mehr Golf zu spielen.

Giuseppe
quelle