Pannaan nyt tämäkin koodi jakoon, vaikka odotettavissa on kateutta ja purnausta.
Jos et usko, että koodi toteuttaa kaikki kunnan säännöt, tee main:iin testikoodi, helppoa, kun nyt main tutkii vain vaihdantalain toteutumista.
/****************************************************************************** * Jouni Aro * * 20.2.20 * * * * The multidimensional algebra * * * * The possibilities of algebra is associated with a sort of tool, with which * * parallel causations can be harnessed and approximated into compact * * functions. When you draw a roughly smooth curve into the screen with the * * mouse, the algebra is able to approximate the given set of dots into a * * short function. The craft of algebra is in its ability to describe * * phenomenons and the nature. * * * * What if you draw curves of different colors into the screen? The white set * * of dots is the temperature, the yellow is the direction of the wind, the * * orange is the time of year, the blue would represent the humidity etc. * * Every set of dots can be described as an individual phenomenon, but what * * about the function, that alone would describe the causations of the whole * * set? * * * * The relation of the function can also be asymmetrical. When the input of * * the weather function can be n number of initations, the output can only be * * a relatively blurry parameter between 0 to 1, whether you should go fishing * * or not. * * * * *** * * * * The computer can be programmed to study mathematics. With the solid axioms * * of algebra and the consistent development of series, the bit monster can * * immerse itself further into geometry and symmetry. In the final solution * * the computer was partly forced to abandon the symbolic concepts which * * include the positive, the negative, the real and the imaginary. The algebra * * doesn't need to share its opinion on in what direction the monitor is put * * on the table. * * * * The algorithms of algebra are constant and its objects can describe * * whatsoever. One of the most beautiful formulae of the algebra is the * * polynomial adaption of the smallest sum of squares: * * * * If the object that is placed in the formula is valid and balanced with * * itself, the algebra will process the substance with rational results, * * otherwise not. A good example is the quaternion. The algebra does not * * editorialise to the inner calculation logic of the quartenion. * * * * *** * * * * i^2 is -1, but what is the object that gives the result of -i when squared? * * There's an interesting connection between the companions of the i, * * |e^(xi+yj+zk+...)| = 1. The multidimensional algebra also connects the pi * * and the number one. pi^z limit on the real axis approaches the value of 1, * * when the z-unit vector is excited further and further from the real axis. * * At the same time, the sum of the absolute values of the axises approaches * * to the pi. * * * * Next up are the simple calculation rules of the recomplex number with the * * C++ examples. The class itself can do almost as such in the MatLab -tool. * * * * ~ recomplex * * ~ order * * ~ specification of class * * ~ addition * * ~ subtraction * * ~ multiplication * * ~ products of vector basis * * ~ division * * ~ abs(recomplex) * * ~ pow(recomplex, int) * * ~ ln(recomplex) * * ~ exp(recomplex) * * ~ pow(recomplex, double) * * ~ pow(recomplex, recomplex) * * ~ sin(recomplex) * * ~ cos(recomplex) * * ~ the port example of the smallest sum of squares * * ~ the connection of the pi and number one * * ~ code attachment * * * * *** * * * * The complex number is a real subset of the recomplex number. Real numbers * * and imaginary numbers for their parts are real subsets of the complex * * number. The prefix "re" points to the feedback, repeating the same basic * * logic in the next scale. * * * * ~~~ order ~~~ * * * * The order goes in power of two, 1, 2, 4, ..., 2^n. The real numbers are * * adjusted with the order of 1. The complex numbers come into effect with the * * order of 2. In the order of 4 the complex numbers get a complex object on * * their side, functioning by the stopper logic. The union of these makes * * simple function maps for the parallel causation phenomena possible, among * * other things. * /*****************************************************************************/ /****************************************************************************** * ~~~ class specification ~~~ * * * * C++ implemented an operator load, when especially the calculation * * operations of the numeric structures on the computer don't differ from the * * regular straightforwardness. * ******************************************************************************/ #include <math.h> #include <stdio.h> #include <conio.h> #include <memory.h> #include <stdlib.h> class recomplex { public: /************************************ * The dimension should be: * * 1, 2, 4, ..., 2^n * * DIMENSION 1 -> real numbers * * DIMENSION 2 -> complex numbers * * DIMENSION 4 -> 4D complex numbers * * DIMENSION 8 -> 8D complex numbers * ************************************/ #define DIMENSION 4 double e[DIMENSION]; recomplex(void); ~recomplex(void); recomplex(double*); recomplex(int); recomplex(int, int); recomplex(int, int, int, int); recomplex(double); recomplex(double, double); recomplex(double, double, double, double); friend void print(recomplex); friend double abs(recomplex); friend recomplex sqrt(recomplex); friend recomplex operator-(recomplex); friend recomplex operator+(recomplex, recomplex); friend recomplex operator-(recomplex, recomplex); friend recomplex operator*(recomplex, recomplex); friend recomplex operator/(recomplex, recomplex); friend recomplex operator*(recomplex, double); friend recomplex operator/(recomplex, double); friend recomplex operator*(double, recomplex); friend recomplex operator/(double, recomplex); friend recomplex pow(recomplex, int); }; /****************************************************************************** * ~~~ addition ~~~ * * * * If you're unfamiliar with the C, a few words about it. The basis of C is * * to save sheets. Every excessive and unnecessary rambling is eliminated * * from the syntax. The operators if C are very effective. For example +=, *=, * * /=, ... perform the operation in question and after that the assignment. * * a.e[n]+=b.e[n] means the same that a.e[n]=a.e[n]+b.e[n]. Because the latter * * style consumes more of paper, the shorter version of the same thing * * a.e[n]+=b.e[n] is used in the C. * ******************************************************************************/ recomplex operator+(recomplex a, recomplex b) { for (int n=0; n<DIMENSION; n++) a.e[n]+=b.e[n]; return a; } /****************************************************************************** * For example: * * * * recomplex A(2, 3, -5, 7); * * recomplex B(7, 5, -3, 2); * * print(A + B); * * * * Prints a value of 9 + 8i - 8j + 9k. * ******************************************************************************/ /****************************************************************************** * ~~~ subtraction ~~~ * ******************************************************************************/ recomplex operator-(recomplex a, recomplex b) { for (int n=0; n<DIMENSION; n++) a.e[n]-=b.e[n]; return a; } /****************************************************************************** * For exsample: * * * * recomplex A(2, 3, 5, 7); * * recomplex B(7, 5, 3, 2); * * print(A-B); * * * * The value of the difference is -5 - 2i + 2j + 5k * ******************************************************************************/ /****************************************************************************** * ~~~ multiplication ~~~ * * * * The multiplication first requires a product definition of the base vectors. * * Looping the product does not editorialise to the formal positive, negative, * * real, imaginary, etc. definitions. The product of the multiplication is * * returned to the sum form so far only in the last loop. * ******************************************************************************/ recomplex operator*(recomplex a, recomplex b) { int *BasisVectors(void); int *R = BasisVectors(); int i, j, n=DIMENSION*2; double x[DIMENSION*2]; double y[DIMENSION*2]; double t[DIMENSION*2]; for (i=j=0; i<DIMENSION; i++, j+=2) x[j+1]=y[j+1]=t[j]=t[j+1]=0.0, x[j]=a.e[i], y[j]=b.e[i]; for (i=0; i<n; i++) for (j=0; j<n; j++) t[R[i*n+j]]+=x[i]*y[j]; for (i=j=0; i<DIMENSION; i++, j+=2) t[i]=(double)(t[j] - t[j+1]); return *(recomplex*) t; } /****************************************************************************** * For exsample: * * * * recomplex A(2, 3, 5, 7); * * recomplex B(7, 5, 3, 2); * * print(A * B); * * * * The outcome is -32 + 32i + 87k. The multiplication is commutative. If A and * * B are complex, the outcome is complex. * * * * recomplex A(2, 3, 0, 0); * * recomplex B(7, 5, 0, 0); * * print(A * B); * * * * You'll get the output -1 + 31i * ******************************************************************************/ /****************************************************************************** * ~~~ the products of the base vectors ~~~ * * * * The products of the base vectors are generated by the recursive function. * * All the products of the base vectors associated with the order of 4 are * * assembled in the chart: * * * * # 1 | -1 | i | -i | j | -j | k | -k | * * ############################################ * * 1 # 1 | -1 | i | -i | j | -j | k | -k | * * ---#----|----|----|----|----|----|----|----| * * -1 # -1 | 1 | -i | i | -j | j | -k | k | * * ---#----|----|----|----|----|----|----|----| * * i # i | -i | -1 | 1 | k | -k | -j | j | * * ---#----|----|----|----|----|----|----|----| * * -i # -i | i | 1 | -1 | -k | k | j | -j | * * ---#----|----|----|----|----|----|----|----| * * j # j | -j | k | -k | i | -i | -1 | 1 | * * ---#----|----|----|----|----|----|----|----| * * -j # -j | j | -k | k | -i | i | 1 | -1 | * * ---#----|----|----|----|----|----|----|----| * * k # k | -k | -j | j | -1 | 1 | -i | i | * * ---#----|----|----|----|----|----|----|----| * * -k # -k | k | j | -j | 1 | -1 | i | -i | * * ---#----|----|----|----|----|----|----|----| * * * * For example ijk = -i, and particularly k^2 = -i. Commutativity remains, * * because still: * * * * ikj = jik = jki = kij = kji = -i * ******************************************************************************/ void GenerateTheProductsOfTheBaseVectors(int *R, int n) { int I=0, J=n-1, i; int X=n, Y=J+n, j; int k=DIMENSION*2; for (i=I; i<=J; i++) for (j=I; j<=J; j++) { R[i*k+X+j]=R[i*k+j]+X; R[(X+j)*k+i]=R[i*k+j]+X; R[(Y-i)*k+Y-j]=J-R[i*k+j]; } if (n+n < DIMENSION*2) { GenerateTheProductsOfTheBaseVectors(R, n+n); } } int* BasisVectors(void) { static int R[DIMENSION*DIMENSION*4]={-1}; if (R[0] == -1) { int FirstThereWasZero=0; R[0x00]=(int)FirstThereWasZero; GenerateTheProductsOfTheBaseVectors(R, 1); } return R; } /****************************************************************************** * ~~~ division ~~~ * * * * In division, a complex conjugate is produced for the divisor, with which * * the divisor and the dividend are multiplied. In the end of every loop step * * there's a half less of axel values of the divisor (the rest are summed up * * to nulls). When the value of the divisor is real, the elements of the * * dividend are divided by the value. For example 2 + 3i + 5j + 7k the mounted * * number is 2 + 3i - 5j - 7k. * ******************************************************************************/ recomplex operator/(recomplex x, recomplex y) { recomplex z; for (int i, n=DIMENSION; n>=2; n/=2) { for (z=y, i=n/2; i<n; i++) z.e[i] = -z.e[i]; x=x*z; y=y*z; } return x/y.e[0]; } /****************************************************************************** * For exsample * * * * recomplex A(2, 3, 5, 7); * * recomplex B(7, 5, 3, 2); * * print(A / B); * * * * The product of division is 0.5488136 - 0.1575801i + 0.7181670j + 0.3977540k.* * If A and B are complex, the result is complex: * * * * recomplex A(2, 3, 0, 0); * * recomplex B(7, 5, 0, 0); * * * print(A / B); * * * * The output of the value is 0.3918919 + 0.1486486i * ******************************************************************************/ recomplex operator/(recomplex x, double k) { for (int i=0; i<DIMENSION; i++) x.e[i]/=(double)k; return x; } recomplex operator*(double k, recomplex x) { for (int i=0; i<DIMENSION; i++) x.e[i]*=(double)k; return x; } recomplex operator*(recomplex x, double k) { return k*x; } recomplex operator/(double k, recomplex x) { return recomplex(k)/x; } /****************************************************************************** * ~~~ abs(recomplex) ~~~ * * * * In the absolute value x is multiplied by its complex conjugate repeatedly. * * The value then returns the breaking power collected by j. * ******************************************************************************/ double abs(recomplex x) { recomplex z; double r, c; int i, j=1, n; for (n=DIMENSION; n>=2; n/=2, j+=j) { for (z=x, i=n/2; i<n; i++) z.e[i]=-z.e[i]; x=x*z; } r=fabs(x.e[0]); c=1.0/(double)j; return pow(r, c); } /****************************************************************************** * For exsample: * * * * recomplex A(2, 3, 5, 7); * * recomplex B(7, 5, 3, 2); * * printf("%25.20f\n", abs(A/B)); * * * * The double of the C is precise, because the absolute value of the result is * * exactly 1.00000000000000000000 * ******************************************************************************/ /****************************************************************************** * ~~~ pow(recomplex, int) ~~~ * * * * The limit of x^0 is 1, when x approaches the null. Calculators have * * variable interpretations, whether the 0^0 is about 1 or is the calculation * * aborted and is the output of the status an error. In the polynomial * * adaption of the smallest sum of squares x^0 has to be 1 also with the x's * * limit null. An exponent n can be a positive or a negative integer. * ******************************************************************************/ recomplex pow(recomplex x, int n) { if (n) { recomplex t=x; int i=n<0? -n: n; for (--i; i; i--) t=t*x; return n>0? t: recomplex(1) / t; } else { return recomplex(1); } } /****************************************************************************** * ~~~ print(recomplex) ~~~ * * * * Somewhere around 20 element, when the marking of the base vector is bigger * * than 'z', the marking of the base vector can be a symbol differentiated * * from the letter. In the calculations themselves the symbols are not * * editorialised. * ******************************************************************************/ void print(recomplex a) { printf("%0.9f ", a.e[0]); for (int n=1; n<DIMENSION; n++) printf("%c %0.9f%c ", a.e[n]<0? '-': '+', fabs(a.e[n]), 'h'+n); printf("\n"); } /****************************************************************************** * These functions still need to be explored * * * * ~~~ ln(recomplex) ~~~ * * ~~~ exp(recomplex) ~~~ * * ~~~ pow(recomplex, double) ~~~ * * ~~~ pow(recomplex, recomplex) ~~~ * * ~~~ sin(recomplex) ~~~ * * ~~~ cos(recomplex) ~~~ * ******************************************************************************/ recomplex::recomplex(void) { memset(this, 0, sizeof(recomplex)); } recomplex::~recomplex(void) { } recomplex::recomplex(double *a) { memset(this, 0, sizeof(recomplex)); for (int i=0x00; i<DIMENSION; i++) e[i]=a[i]; } recomplex::recomplex(int a) { memset(this, 0, sizeof(recomplex)); e[0]=(double)a; } recomplex::recomplex(double a) { memset(this, 0, sizeof(recomplex)); e[0]=(double)a; } recomplex::recomplex(int a, int b) { memset(this, 0, sizeof(recomplex)); e[0]=(double)a; e[1]=(double)b; } recomplex::recomplex(double a, double b) { memset(this, 0, sizeof(recomplex)); e[0]=(double)a; e[1]=(double)b; } recomplex::recomplex(int a, int b, int c, int d) { memset(this, 0, sizeof(recomplex)); e[0]=(double)a; e[1]=(double)b; e[2]=(double)c; e[3]=(double)d; } recomplex::recomplex(double a, double b, double c, double d) { memset(this, 0, sizeof(recomplex)); e[0]=(double)a; e[1]=(double)b; e[2]=(double)c; e[3]=(double)d; } /*************************************************** * Returns the square root of parameter c. Function * * has been implemented in Newton's processing. * ***************************************************/ recomplex sqrt(recomplex c) { int i; recomplex x; double R[DIMENSION]; for (i=0; i<DIMENSION; i++) { int sgn = rand()%0x02? -1: +1; R[i]=sgn*rand()/(double)32767; } x=recomplex((double*)R); for (i=0; i<1024; i++) x=x-(x*x-c)/(2.0*x); return (recomplex)x; } /********************************************** * the mainfunction tests the income swap law. * **********************************************/ int main(void) { recomplex x(0.12, -0.34, 0.56, -0.78); recomplex y(0.98, -0.76, 0.54, -0.32); printf("x = "); print(x); printf("y = "); print(y); printf("\n"); printf("x*y = "); print(x*y); printf("y*x = "); print(y*x); printf("\n"); getch(); return 0; }
Osaisitko selittää, mikä on tämän järjestelmän ero kompleksilukuihin ja kvaternioihin? Tekninen ero on tietenkin tuo erilainen kantavektorien tulotaulukko, mutta minkä todellisen ongelman tämä järjestelmä ratkaisee kompleksilukuihin ja kvaternioihin verrattuna? Kunnes tähän saadaan selvä vastaus, ”purnaus” on täysin aiheellista ja kateuteen ei ole pienintäkään aihetta. Laajasti tutkitun ja tunnetun järjestelmän käyttö matematiikassa on varmasti parempi vaihtoehto kuin omatekoinen viritelmä, jos ei ole tarjota konkreettista vastaesimerkkiä.
Vuonna 2008 on jo täällä keskusteltu tekemästäsi real2d-viritelmästä, joka oli kommentoijien mukaan matemaattisesti rikkinäinen versio kompleksiluvuista. (Itse et tietysti myöntänyt vikoja silloinkaan.) Nyt näyttää, että olet laajentanut samanlaisen rikkinäisyyden useampaan ulottuvuuteen. Jos olisit tehnyt oikeasti läpimurron matematiikassa, voisit julkaista sen varmasti paremmallakin paikalla kuin Ohjelmointiputkan koodinurkassa.
Tämän kyseisen koodisi laadusta C++:n näkökulmasta on jo aiemmin keskusteltu, ja niitäkään vikoja et ole korjannut, vaan ylläpidät samaa purkkaa.
Vai niin. Fasebookkiin (https://www.facebook.com/jouni.aro.10) tein esimerkin, miten aikajana korreloi 8 käyrän kanssa. Käyrät suodatin pienimmän neliösumman polynomiapproksimaatiolla.
Minulla on perseessä sellainen kutina, että recomplex-lukuja voi käyttää tekoälysovelluksessa.
Solu voisi sisältää 8 käyrää, jotka yhdessä lähettävät yhden reaaliarvon synapsia pitkin seuraavalle solulle, joka taas saa 8 inputtia ja lähettää niiden avulla jälleen yhden herätteen johonkin soluun.
Aika on suhteellista. recomplex-lukujen toteuttamiseen ja tutkimukseen kului likmain 35 vuotta, joten nyt tekoälyohjelman duunaamiseen voi kulua seuraavat 35 vuotta.
Olen sitkeä sissi, enkä anna periksi, vaikka epäonnistuisin miljoona kertaa projektin lomassa. Toki teen kokoajan muitakin projekteja.
Ja vielä kysymykseesi: IO voi olla epäsymmetrinen, esimerkiksi 4 epälineaarista käyrää korreloi 16 käyrän kanssa.
Ja mikä koodissa on mukamas purkkaviritelmää. Minusta se on erittäin kaunis ja korrekti olio.
Kukaan ei tiedä tuosta selityksestäsi, mitä tarkoittaa ”aikajana korreloi 8 käyrän kanssa” tai miten tämä liittyy sinun koodiisi. Pystytkö määrittelemään tämän tavoitteen niin selvästi ja yksiselitteisesti oikealla matematiikan kielellä, että muutkin voisivat ratkaista tehtävän? (Oletko edes määritellyt ongelmaa matemaattisesti, vai heitätkö vain pikseleitä ruudulle ja toteat, että kylläpä nyt ”korreloi käyrän kanssa” hienosti?) Vasta sitten voitaisiin verrata, onko sinun koodistasi jotain hyötyä muihin verrattuna. Toistaiseksi ollaan vielä tilanteessa, jossa heität hienoja termejä ilman mitään konkreettista sisältöä ja kerrot, että jossain toisessa koodissa olet recomplex-lukujesi avulla ratkaissut tämän ongelman hienosti. Tai jos edelleen puhut siitä pienimmän neliösumman menetelmästä tilastotieteessä, niin edellisen keskustelun argumentit pätevät edelleen.
Ylipäänsä ei ole uskottavaa esittää keksintöä matematiikasta, kun perustermitkin ovat väärin. Koodissasi kerrot, että testataan ”income swap law” ([taloudellisen] tulon vaihtamisen laki, mikä lieneekään), kun pitäisi varmaan olla ”commutativity of multiplication” (kertolaskun vaihdannaisuus).
jone2712 kirjoitti:
Ja vielä kysymykseesi: IO voi olla epäsymmetrinen, esimerkiksi 4 epälineaarista käyrää korreloi 16 käyrän kanssa.
En tiedä, mihin kysymykseen tämä nyt mielestäsi vastasi tai mitä edes yrität sanoa. Ehkä tarkoitit, että voidaan ottaa 4 käyrää ja muodostaa 16 yksinkertaisempaa käyrää, jotka approksimoivat niitä tietyllä lukuvälillä. Tämä onnistuu varmasti myös ilman sinun lukujärjestelmääsi.
jone2712 kirjoitti:
Ja mikä koodissa on mukamas purkkaviritelmää. Minusta se on erittäin kaunis ja korrekti olio.
Alkajaisiksi esim. turhat tyypinmuunnokset joka välissä, turha memsetin käyttö, turhat erilliset muodostimet ja turha aputaulukoiden käyttö. Näitä luettelin ja korjasin jo viimeksi. Vielä parempiakin ratkaisuja on, mutta turha niitä on nyt tässä selittää, kun et taida olla kuitenkaan kiinnostunut.
Tämä foorumi on niin köyhä, ettei sisällä kuvien julkaisuja. Yksi kuva kertoo yleensä enemmän kuin 1000 sanaa. Joten kirjoitelkaa vain toisillenne eksakteilla koodilla ja merkinnöillä, niin kaikki voivat olla iloisia.
Kun moitit koodia, että siinä on turhia elementtejä, koodaa näytöksi yksikin funktion pätkä, joka on merkittävästi tehokkaampi, ja nopeampi. sqrt2 voisi optimoida, mutta ei se ole merkittävää vielä tässä vaiheessa, kun koodi on vasta demo.
Ja aloituspuheenvuorossa oli sanottu vastaus kysymykseesi. Minun mielestäni purnauksessasi on katkeruutta mukana. Sanoin silloin ja silloin, että pitää olla niin ja näin, vattu oliossa on yli 30 vuoden kokemus. Sinäkin olet varmaan tehnyt shakkialgoritmin 64 assemblerilla?
Mistä voi ladata Julia kielen
T: nero
jone2712 kirjoitti:
Yksi kuva kertoo yleensä enemmän kuin 1000 sanaa.
Voisit kokeilla kuitenkin niitä 1000 sanaa. Jos et osaa edes selittää asiaa, et voi keksiä uudenlaista toimivaa matematiikkaa. Kuvapalveluksi käy vaikka https://snipboard.io/, mutta kuva ei ole mikään matemaattisen ongelman tai keksinnön kuvaus. Et voi tosissasi väittää, että uusi algebran muoto toimii, koska tämä kuva näyttää hyvältä.
jone2712 kirjoitti:
purnauksessasi on katkeruutta – –
yli 30 vuoden kokemus – –
shakkialgoritmin 64 assemblerilla
Oho, vuosi sitten oli 27 vuoden kokemus ja nyt jo yli 30! Kokemus ei auta, jos ei kehity. Vuoden takaiseen nähden koodisi ei ole kehittynyt, joten ei taida olla kertynyt myöskään hyödyllistä kokemusta asiasta. Argumentitkin näköjään menevät vanhalle linjalle, eli kun et osaa perustella, alat selittää purnauksesta ja katkeruudesta ja leuhkia jollain ikivanhoilla shakkikoodeilla. Jos ”keskustelu” jatkuu sillä linjalla, kannattaa varmaan lukita koko ketju.
jone2712 kirjoitti:
Kun moitit koodia, että siinä on turhia elementtejä, koodaa näytöksi yksikin funktion pätkä, joka on merkittävästi tehokkaampi, ja nopeampi.
Edellisessä keskustelussa näytin esimerkiksi, miten kertolaskusi onnistuu paljon vähemmillä toimituksilla ja miten kahdeksan erillistä muodostinta korvataan yhdellä. Sitä paitsi en kritisoinut koodin tehokkuutta vaan muuta laatua.
No, oletetaan että saat tehtävän:
Optimoi 8 polynomifunktioon aikajana 0, 1, 2, …, 31. Polynomifunktioissa on kohinaa, joten joudut käyttämään pienimmän neliösumman polynomiapproksimaatiota seuraavasti:
polynomi(*poly1, *x1, *y1) polynomi(*poly2, *x2, *y2) polynomi(*poly3, *x3, *y3) polynomi(*poly4, *x4, *y4) polynomi(*poly5, *x5, *y5) polynomi(*poly6, *x6, *y6) polynomi(*poly7, *x7, *y7) polynomi(*poly8, *x8, *y8)
x on aika, ja y on sitä vastaavan polynomin arvot.
Siis nyt sinulla on toisessa funktiossa esimerkiksi piirto:
double x1=poly1(x1) double x2=poly2(x2) double x3=poly3(x3) double x4=poly4(x4) double x5=poly5(x5) double x6=poly6(x6) double x7=poly7(x7) double x8=poly8(x8) plot(1, y1) plot(2, y2) plot(3, y3) plot(4, y4) plot(5, y5) plot(6, y6) plot(7, y7) plot(8, y8), jne.
Pointti on kuitenkin se, että joudut käpristelemään 8 polynomifunktion kanssa. Tulee helposti painovirhe.
16D-recomplex luvulla homma on vähän yksinkertaisempi:
mk_poly(poly, in, out)
Nyt recomplex polynomi pitää sisällään koko funktioparven.
Piirto funktio kävisi luupissa jokseenkin seuraavasti:
void piirra(recomplex *poly) { recomplex y=poly.fx(x++) plot(x, y) }
Mutta tästä on seikkaperäisempi esimerkki minun fasebookissa osoitteessa:
https://www.facebook.com/groups/6784853553/
Ja miksi 16D-recomplex polynomi?
Outputin arvot on sijoitettava parillisille akseleille. Kompleksiset parittomat akselit saavat sitten arvonsa polynomin approksimoinnissa.
oletko ajatellut mitä teet rahoilla, kun saanet tästä tänä vuonna noopelin palkinnon?
Matematiikassa ei jaeta nobeleita.
kai tästä jotain fyrkkaa on parempi saada jos tätä olet 30 vuotta vääntänyt
jone2712 kirjoitti:
Optimoi 8 polynomifunktioon aikajana 0, 1, 2, …, 31. – – Pointti on kuitenkin se, että joudut käpristelemään 8 polynomifunktion kanssa. Tulee helposti painovirhe.
Tiivistän vielä, että puhuttaisiin varmasti samasta asiasta.
Esittämäsi tehtävä on siis ilmeisesti se, että on annettu jokin kahdeksan sarakkeen (plus indeksisarake) datajoukko, jossa siis jokaisesta sarakkeesta piirretään erillinen viiva kuvaan ja jokaista saraketta ”optimoidaan polynomifunktioon” eli yritetään kuvata polynomilla. Pienimmän neliösumman menetelmä on yksi tapa polynomin kertoimien etsimiseen.
Esitit siihen oikean ratkaisun, eli jokaiselle sarakkeelle tehdään oma polynomi. Mielestäsi ratkaisu on huono, koska ”helposti tulee painovirhe”. No tämä puute ratkeaa, kun opettelee käyttämään vaikka taulukoita tai valarray-tyyppiä eikä copypastea joka sarakkeelle. Muita vikoja et valittanut, joten ilmeisesti matematiikan puolesta ratkaisu kelpaa. Eli myönnät itsekin, että ratkaisuun ei tarvita recomplex-lukuja.
jone2712 kirjoitti:
Polynomifunktioissa on kohinaa,
Polynomifunktioissa ei ole kohinaa. Tarkoitat ehkä, että annetussa datassa on ”kohinaa” verrattuna siitä laskettavaan polynomiin.
jone2712 kirjoitti:
Mutta tästä on seikkaperäisempi esimerkki minun fasebookissa
Linkittämäsi postaus sisälsi kuvan käyristä ja noin 50 sanaa, joten ei se kovin seikkaperäinen ole. Kirjoitat Facebook-postauksessa: ”Vielä kiinnostaa, minkälaiset polynomin kertoimet ovat.” En ymmärrä, miten voit selittää huippuhienosta ratkaisustasi, jos et saa siitä selville edes polynomin kertoimia.
jone2712 kirjoitti:
16D-recomplex luvulla homma on vähän yksinkertaisempi: ... Nyt recomplex polynomi pitää sisällään koko funktioparven.
Jos polynomi on muotoa axⁿ + ... + bx² + cx + d, missä x on tavallinen reaaliluku (kuten kasvava indeksi) ja kertoimet ovat yksittäisiä recomplex-lukuja, olet sinänsä oikeassa siinä, että tässä saadaan myös tulokseksi recomplex-luku, jossa on usean polynomin ratkaisut. Harmi (tai onni), että reaaliluvun kanssa kertolaskussa recomplex-luku ei sisällä mitään uutta moniulotteista algebraa vaan on vain valtava purkkapallo tavallisen lukutaulukon ympärillä.
Jos taas polynomissasi on jossain kohti kahden (useampia ulottuvuuksia sisältävän) recomplex-luvun kertolasku, tilanne on aivan toinen ja pitäisi nähdä se todellinen käyttämäsi polynomikoodi, jotta voisi sanoa, mitä siinä matemaattisesti tapahtuu. Kuitenkin lähtökohtaisesti kertolaskussa sinulla dimensiot muuttuvat toisikseen ja tuloksena on todennäköisesti jotain matemaattisesti järjetöntä, kun kuitenkin tarkoitus oli, että jokaisessa recomplex-luvun ”dimensiossa” käsitellään yhtä syötetaulukon saraketta riippumatta muista sarakkeista.
jone2712 kirjoitti:
Kun moitit koodia, että siinä on turhia elementtejä, koodaa näytöksi yksikin funktion pätkä, joka on merkittävästi tehokkaampi, ja nopeampi.
Tässä on sinulle nykyaikainen (C++20) muodostin eli konstruktori, joka suoraan korvaa kaikki paitsi pointteriversion. Säästät tällä noin 60 riviä koodia.
template <typename... T> recomplex(T... init_data): e {init_data...} {}
Viime vuoden keskustelusta löytyy versio, joka toimii hieman eri syntaksilla vanhemmilla C++-versioilla.
Ota huomioon, että koodi kommentteineen on tarkoitettu tavallisille matematiikan osaajille, joilla on vähän koulutusta C-ohjelmointiin. Tuota koodia pystyy lukemaan vähemmälläkin C-kokemuksella. Ei mitään tempalete virityksiä. Eikä se funktiokäyräparveen optimointi ihan noin mennyt. Tarkennan sitä, kun on aikaa.
Mutta käyräparven laskennassa sekä parilliset että parittomat akselit asettuvat. Ei siis niin, että vain parilliset akselit asettuisivat.
Tuossa on input ja output taulukot. Tee niille pienimmän neliösumman polynomiapproksimaatio, ja tarkastele graafisesti, kuinka hyvin approksimaatio toimii.
recomplex input[68]= { recomplex(0.0000), recomplex(0.0147), recomplex(0.0294), recomplex(0.0441), recomplex(0.0588), recomplex(0.0735), recomplex(0.0882), recomplex(0.1029), recomplex(0.1176), recomplex(0.1324), recomplex(0.1471), recomplex(0.1618), recomplex(0.1765), recomplex(0.1912), recomplex(0.2059), recomplex(0.2206), recomplex(0.2353), recomplex(0.2500), recomplex(0.2647), recomplex(0.2794), recomplex(0.2941), recomplex(0.3088), recomplex(0.3235), recomplex(0.3382), recomplex(0.3529), recomplex(0.3676), recomplex(0.3824), recomplex(0.3971), recomplex(0.4118), recomplex(0.4265), recomplex(0.4412), recomplex(0.4559), recomplex(0.4706), recomplex(0.4853), recomplex(0.5000), recomplex(0.5147), recomplex(0.5294), recomplex(0.5441), recomplex(0.5588), recomplex(0.5735), recomplex(0.5882), recomplex(0.6029), recomplex(0.6176), recomplex(0.6324), recomplex(0.6471), recomplex(0.6618), recomplex(0.6765), recomplex(0.6912), recomplex(0.7059), recomplex(0.7206), recomplex(0.7353), recomplex(0.7500), recomplex(0.7647), recomplex(0.7794), recomplex(0.7941), recomplex(0.8088), recomplex(0.8235), recomplex(0.8382), recomplex(0.8529), recomplex(0.8676), recomplex(0.8824), recomplex(0.8971), recomplex(0.9118), recomplex(0.9265), recomplex(0.9412), recomplex(0.9559), recomplex(0.9706), recomplex(0.9853), }; recomplex output[68]= { recomplex(+0.2323, 0.0, +0.2634, 0.0, -0.3794, 0.0, +0.3356, 0.0, -0.2773, 0.0, +0.3614, 0.0, -0.2055, 0.0, 0.0, 0.0), recomplex(+0.2093, 0.0, +0.2297, 0.0, -0.3078, 0.0, +0.3020, 0.0, -0.2945, 0.0, +0.3587, 0.0, -0.0906, 0.0, 0.0, 0.0), recomplex(+0.1869, 0.0, +0.1989, 0.0, -0.2443, 0.0, +0.2747, 0.0, -0.2892, 0.0, +0.3461, 0.0, -0.0239, 0.0, 0.0, 0.0), recomplex(+0.1651, 0.0, +0.1707, 0.0, -0.1884, 0.0, +0.2529, 0.0, -0.2663, 0.0, +0.3271, 0.0, +0.0086, 0.0, 0.0, 0.0), recomplex(+0.1440, 0.0, +0.1452, 0.0, -0.1395, 0.0, +0.2359, 0.0, -0.2302, 0.0, +0.3043, 0.0, +0.0181, 0.0, 0.0, 0.0), recomplex(+0.1235, 0.0, +0.1223, 0.0, -0.0970, 0.0, +0.2232, 0.0, -0.1846, 0.0, +0.2799, 0.0, +0.0135, 0.0, 0.0, 0.0), recomplex(+0.1037, 0.0, +0.1019, 0.0, -0.0606, 0.0, +0.2140, 0.0, -0.1328, 0.0, +0.2557, 0.0, +0.0011, 0.0, 0.0, 0.0), recomplex(+0.0845, 0.0, +0.0838, 0.0, -0.0297, 0.0, +0.2080, 0.0, -0.0775, 0.0, +0.2329, 0.0, -0.0140, 0.0, 0.0, 0.0), recomplex(+0.0659, 0.0, +0.0680, 0.0, -0.0039, 0.0, +0.2044, 0.0, -0.0213, 0.0, +0.2125, 0.0, -0.0285, 0.0, 0.0, 0.0), recomplex(+0.0480, 0.0, +0.0544, 0.0, +0.0173, 0.0, +0.2030, 0.0, +0.0341, 0.0, +0.1950, 0.0, -0.0403, 0.0, 0.0, 0.0), recomplex(+0.0307, 0.0, +0.0429, 0.0, +0.0344, 0.0, +0.2031, 0.0, +0.0869, 0.0, +0.1807, 0.0, -0.0480, 0.0, 0.0, 0.0), recomplex(+0.0140, 0.0, +0.0334, 0.0, +0.0477, 0.0, +0.2046, 0.0, +0.1360, 0.0, +0.1697, 0.0, -0.0510, 0.0, 0.0, 0.0), recomplex(-0.0020, 0.0, +0.0258, 0.0, +0.0576, 0.0, +0.2070, 0.0, +0.1803, 0.0, +0.1621, 0.0, -0.0494, 0.0, 0.0, 0.0), recomplex(-0.0174, 0.0, +0.0200, 0.0, +0.0646, 0.0, +0.2100, 0.0, +0.2191, 0.0, +0.1576, 0.0, -0.0436, 0.0, 0.0, 0.0), recomplex(-0.0321, 0.0, +0.0159, 0.0, +0.0690, 0.0, +0.2133, 0.0, +0.2520, 0.0, +0.1558, 0.0, -0.0341, 0.0, 0.0, 0.0), recomplex(-0.0462, 0.0, +0.0135, 0.0, +0.0711, 0.0, +0.2167, 0.0, +0.2787, 0.0, +0.1565, 0.0, -0.0218, 0.0, 0.0, 0.0), recomplex(-0.0597, 0.0, +0.0127, 0.0, +0.0713, 0.0, +0.2200, 0.0, +0.2993, 0.0, +0.1590, 0.0, -0.0075, 0.0, 0.0, 0.0), recomplex(-0.0725, 0.0, +0.0132, 0.0, +0.0699, 0.0, +0.2230, 0.0, +0.3138, 0.0, +0.1631, 0.0, +0.0079, 0.0, 0.0, 0.0), recomplex(-0.0847, 0.0, +0.0152, 0.0, +0.0672, 0.0, +0.2255, 0.0, +0.3224, 0.0, +0.1682, 0.0, +0.0235, 0.0, 0.0, 0.0), recomplex(-0.0962, 0.0, +0.0184, 0.0, +0.0635, 0.0, +0.2275, 0.0, +0.3256, 0.0, +0.1739, 0.0, +0.0388, 0.0, 0.0, 0.0), recomplex(-0.1071, 0.0, +0.0228, 0.0, +0.0590, 0.0, +0.2287, 0.0, +0.3237, 0.0, +0.1797, 0.0, +0.0530, 0.0, 0.0, 0.0), recomplex(-0.1174, 0.0, +0.0283, 0.0, +0.0541, 0.0, +0.2292, 0.0, +0.3175, 0.0, +0.1852, 0.0, +0.0657, 0.0, 0.0, 0.0), recomplex(-0.1270, 0.0, +0.0347, 0.0, +0.0488, 0.0, +0.2289, 0.0, +0.3073, 0.0, +0.1902, 0.0, +0.0766, 0.0, 0.0, 0.0), recomplex(-0.1360, 0.0, +0.0421, 0.0, +0.0436, 0.0, +0.2277, 0.0, +0.2940, 0.0, +0.1942, 0.0, +0.0856, 0.0, 0.0, 0.0), recomplex(-0.1444, 0.0, +0.0502, 0.0, +0.0384, 0.0, +0.2256, 0.0, +0.2782, 0.0, +0.1970, 0.0, +0.0925, 0.0, 0.0, 0.0), recomplex(-0.1521, 0.0, +0.0591, 0.0, +0.0336, 0.0, +0.2227, 0.0, +0.2605, 0.0, +0.1985, 0.0, +0.0974, 0.0, 0.0, 0.0), recomplex(-0.1591, 0.0, +0.0686, 0.0, +0.0293, 0.0, +0.2189, 0.0, +0.2415, 0.0, +0.1986, 0.0, +0.1005, 0.0, 0.0, 0.0), recomplex(-0.1656, 0.0, +0.0786, 0.0, +0.0257, 0.0, +0.2143, 0.0, +0.2221, 0.0, +0.1971, 0.0, +0.1021, 0.0, 0.0, 0.0), recomplex(-0.1714, 0.0, +0.0891, 0.0, +0.0229, 0.0, +0.2090, 0.0, +0.2027, 0.0, +0.1941, 0.0, +0.1024, 0.0, 0.0, 0.0), recomplex(-0.1765, 0.0, +0.0999, 0.0, +0.0209, 0.0, +0.2029, 0.0, +0.1840, 0.0, +0.1896, 0.0, +0.1019, 0.0, 0.0, 0.0), recomplex(-0.1811, 0.0, +0.1110, 0.0, +0.0200, 0.0, +0.1962, 0.0, +0.1664, 0.0, +0.1838, 0.0, +0.1009, 0.0, 0.0, 0.0), recomplex(-0.1849, 0.0, +0.1222, 0.0, +0.0202, 0.0, +0.1890, 0.0, +0.1505, 0.0, +0.1767, 0.0, +0.0998, 0.0, 0.0, 0.0), recomplex(-0.1882, 0.0, +0.1335, 0.0, +0.0215, 0.0, +0.1812, 0.0, +0.1366, 0.0, +0.1685, 0.0, +0.0990, 0.0, 0.0, 0.0), recomplex(-0.1908, 0.0, +0.1447, 0.0, +0.0240, 0.0, +0.1731, 0.0, +0.1251, 0.0, +0.1594, 0.0, +0.0987, 0.0, 0.0, 0.0), recomplex(-0.1928, 0.0, +0.1559, 0.0, +0.0278, 0.0, +0.1647, 0.0, +0.1163, 0.0, +0.1497, 0.0, +0.0993, 0.0, 0.0, 0.0), recomplex(-0.1941, 0.0, +0.1668, 0.0, +0.0329, 0.0, +0.1562, 0.0, +0.1102, 0.0, +0.1396, 0.0, +0.1010, 0.0, 0.0, 0.0), recomplex(-0.1948, 0.0, +0.1774, 0.0, +0.0392, 0.0, +0.1475, 0.0, +0.1071, 0.0, +0.1292, 0.0, +0.1040, 0.0, 0.0, 0.0), recomplex(-0.1948, 0.0, +0.1877, 0.0, +0.0467, 0.0, +0.1388, 0.0, +0.1069, 0.0, +0.1190, 0.0, +0.1082, 0.0, 0.0, 0.0), recomplex(-0.1942, 0.0, +0.1974, 0.0, +0.0555, 0.0, +0.1303, 0.0, +0.1095, 0.0, +0.1090, 0.0, +0.1138, 0.0, 0.0, 0.0), recomplex(-0.1930, 0.0, +0.2066, 0.0, +0.0655, 0.0, +0.1219, 0.0, +0.1149, 0.0, +0.0995, 0.0, +0.1206, 0.0, 0.0, 0.0), recomplex(-0.1911, 0.0, +0.2151, 0.0, +0.0765, 0.0, +0.1138, 0.0, +0.1227, 0.0, +0.0906, 0.0, +0.1286, 0.0, 0.0, 0.0), recomplex(-0.1886, 0.0, +0.2228, 0.0, +0.0886, 0.0, +0.1061, 0.0, +0.1328, 0.0, +0.0825, 0.0, +0.1374, 0.0, 0.0, 0.0), recomplex(-0.1855, 0.0, +0.2297, 0.0, +0.1016, 0.0, +0.0988, 0.0, +0.1447, 0.0, +0.0752, 0.0, +0.1468, 0.0, 0.0, 0.0), recomplex(-0.1817, 0.0, +0.2356, 0.0, +0.1154, 0.0, +0.0920, 0.0, +0.1581, 0.0, +0.0688, 0.0, +0.1567, 0.0, 0.0, 0.0), recomplex(-0.1773, 0.0, +0.2405, 0.0, +0.1299, 0.0, +0.0858, 0.0, +0.1724, 0.0, +0.0634, 0.0, +0.1666, 0.0, 0.0, 0.0), recomplex(-0.1722, 0.0, +0.2443, 0.0, +0.1449, 0.0, +0.0802, 0.0, +0.1872, 0.0, +0.0587, 0.0, +0.1763, 0.0, 0.0, 0.0), recomplex(-0.1665, 0.0, +0.2468, 0.0, +0.1602, 0.0, +0.0752, 0.0, +0.2019, 0.0, +0.0548, 0.0, +0.1854, 0.0, 0.0, 0.0), recomplex(-0.1602, 0.0, +0.2480, 0.0, +0.1758, 0.0, +0.0708, 0.0, +0.2159, 0.0, +0.0515, 0.0, +0.1939, 0.0, 0.0, 0.0), recomplex(-0.1532, 0.0, +0.2479, 0.0, +0.1913, 0.0, +0.0671, 0.0, +0.2287, 0.0, +0.0485, 0.0, +0.2014, 0.0, 0.0, 0.0), recomplex(-0.1456, 0.0, +0.2462, 0.0, +0.2066, 0.0, +0.0639, 0.0, +0.2397, 0.0, +0.0455, 0.0, +0.2078, 0.0, 0.0, 0.0), recomplex(-0.1374, 0.0, +0.2429, 0.0, +0.2214, 0.0, +0.0613, 0.0, +0.2483, 0.0, +0.0423, 0.0, +0.2132, 0.0, 0.0, 0.0), recomplex(-0.1285, 0.0, +0.2379, 0.0, +0.2356, 0.0, +0.0592, 0.0, +0.2540, 0.0, +0.0385, 0.0, +0.2176, 0.0, 0.0, 0.0), recomplex(-0.1189, 0.0, +0.2312, 0.0, +0.2487, 0.0, +0.0574, 0.0, +0.2565, 0.0, +0.0337, 0.0, +0.2212, 0.0, 0.0, 0.0), recomplex(-0.1088, 0.0, +0.2225, 0.0, +0.2606, 0.0, +0.0560, 0.0, +0.2551, 0.0, +0.0276, 0.0, +0.2243, 0.0, 0.0, 0.0), recomplex(-0.0980, 0.0, +0.2120, 0.0, +0.2709, 0.0, +0.0547, 0.0, +0.2498, 0.0, +0.0198, 0.0, +0.2272, 0.0, 0.0, 0.0), recomplex(-0.0865, 0.0, +0.1993, 0.0, +0.2794, 0.0, +0.0533, 0.0, +0.2402, 0.0, +0.0100, 0.0, +0.2304, 0.0, 0.0, 0.0), recomplex(-0.0744, 0.0, +0.1845, 0.0, +0.2856, 0.0, +0.0517, 0.0, +0.2264, 0.0, -0.0022, 0.0, +0.2345, 0.0, 0.0, 0.0), recomplex(-0.0617, 0.0, +0.1675, 0.0, +0.2893, 0.0, +0.0497, 0.0, +0.2085, 0.0, -0.0169, 0.0, +0.2399, 0.0, 0.0, 0.0), recomplex(-0.0484, 0.0, +0.1482, 0.0, +0.2901, 0.0, +0.0469, 0.0, +0.1867, 0.0, -0.0342, 0.0, +0.2471, 0.0, 0.0, 0.0), recomplex(-0.0344, 0.0, +0.1264, 0.0, +0.2876, 0.0, +0.0430, 0.0, +0.1615, 0.0, -0.0542, 0.0, +0.2564, 0.0, 0.0, 0.0), recomplex(-0.0197, 0.0, +0.1021, 0.0, +0.2814, 0.0, +0.0379, 0.0, +0.1337, 0.0, -0.0765, 0.0, +0.2682, 0.0, 0.0, 0.0), recomplex(-0.0045, 0.0, +0.0752, 0.0, +0.2710, 0.0, +0.0309, 0.0, +0.1043, 0.0, -0.1009, 0.0, +0.2821, 0.0, 0.0, 0.0), recomplex(+0.0115, 0.0, +0.0456, 0.0, +0.2562, 0.0, +0.0219, 0.0, +0.0746, 0.0, -0.1265, 0.0, +0.2976, 0.0, 0.0, 0.0), recomplex(+0.0280, 0.0, +0.0132, 0.0, +0.2363, 0.0, +0.0103, 0.0, +0.0462, 0.0, -0.1525, 0.0, +0.3136, 0.0, 0.0, 0.0), recomplex(+0.0452, 0.0, -0.0221, 0.0, +0.2110, 0.0, -0.0044, 0.0, +0.0212, 0.0, -0.1774, 0.0, +0.3281, 0.0, 0.0, 0.0), recomplex(+0.0630, 0.0, -0.0604, 0.0, +0.1798, 0.0, -0.0228, 0.0, +0.0018, 0.0, -0.1994, 0.0, +0.3382, 0.0, 0.0, 0.0), recomplex(+0.0815, 0.0, -0.1017, 0.0, +0.1421, 0.0, -0.0453, 0.0, -0.0091, 0.0, -0.2164, 0.0, +0.3398, 0.0, 0.0, 0.0), recomplex(+0.1006, 0.0, -0.1462, 0.0, +0.0975, 0.0, -0.0726, 0.0, -0.0083, 0.0, -0.2255, 0.0, +0.3273, 0.0, 0.0, 0.0), };
jone2712 kirjoitti:
Ota huomioon, että koodi kommentteineen on tarkoitettu tavallisille matematiikan osaajille, joilla on vähän koulutusta C-ohjelmointiin. Tuota koodia pystyy lukemaan vähemmälläkin C-kokemuksella...
Hieman päälle vuosikymmenen see-kokemuksella, en kyllä millään muista nähneeni missään mainintaa tuosta class avainsanasta, onko se kenties makro structille jossain?. Tuo public-labeli tuntuu olevan myöskin ihmeellisessä paikassa enkä ole ollenkaan varma onko sen käyttö näin sallittua. friend tuntuu korvaan myös oudolta, onko kenties gcc-laajennus tai muuta sellaista? nuo operator-nimiset funktiot ihmetyttävät myöskin, koska C:n funkkareiden nimissä ei saa esiintyä +-*/ jne, ja tietenkin niiden pitäisi olla pointtereita noin struktin jäseninä. Suosittelen korjaamaan soodisi.
Esimerkkidata kuvaa näköjään suoraan niitä polynomeja ilman mitään ”kohinaa” (paitsi pyöristys neljään desimaaliin), joten tuosta saisi polynomit pulautettua muillakin menetelmillä.
Jotain kertoo myös se, että tuo esimerkkidata on tungettu kahdeksi recomplex-taulukoksi kaikkien turhien nollien kanssa. Tässä on yllä oleva data CSV-muodossa (data.csv).
0,0.2323,0.2634,-0.3794,0.3356,-0.2773,0.3614,-0.2055 0.0147,0.2093,0.2297,-0.3078,0.302,-0.2945,0.3587,-0.0906 0.0294,0.1869,0.1989,-0.2443,0.2747,-0.2892,0.3461,-0.0239 0.0441,0.1651,0.1707,-0.1884,0.2529,-0.2663,0.3271,0.0086 0.0588,0.144,0.1452,-0.1395,0.2359,-0.2302,0.3043,0.0181 0.0735,0.1235,0.1223,-0.097,0.2232,-0.1846,0.2799,0.0135 0.0882,0.1037,0.1019,-0.0606,0.214,-0.1328,0.2557,0.0011 0.1029,0.0845,0.0838,-0.0297,0.208,-0.0775,0.2329,-0.014 0.1176,0.0659,0.068,-0.0039,0.2044,-0.0213,0.2125,-0.0285 0.1324,0.048,0.0544,0.0173,0.203,0.0341,0.195,-0.0403 0.1471,0.0307,0.0429,0.0344,0.2031,0.0869,0.1807,-0.048 0.1618,0.014,0.0334,0.0477,0.2046,0.136,0.1697,-0.051 0.1765,-0.002,0.0258,0.0576,0.207,0.1803,0.1621,-0.0494 0.1912,-0.0174,0.02,0.0646,0.21,0.2191,0.1576,-0.0436 0.2059,-0.0321,0.0159,0.069,0.2133,0.252,0.1558,-0.0341 0.2206,-0.0462,0.0135,0.0711,0.2167,0.2787,0.1565,-0.0218 0.2353,-0.0597,0.0127,0.0713,0.22,0.2993,0.159,-0.0075 0.25,-0.0725,0.0132,0.0699,0.223,0.3138,0.1631,0.0079 0.2647,-0.0847,0.0152,0.0672,0.2255,0.3224,0.1682,0.0235 0.2794,-0.0962,0.0184,0.0635,0.2275,0.3256,0.1739,0.0388 0.2941,-0.1071,0.0228,0.059,0.2287,0.3237,0.1797,0.053 0.3088,-0.1174,0.0283,0.0541,0.2292,0.3175,0.1852,0.0657 0.3235,-0.127,0.0347,0.0488,0.2289,0.3073,0.1902,0.0766 0.3382,-0.136,0.0421,0.0436,0.2277,0.294,0.1942,0.0856 0.3529,-0.1444,0.0502,0.0384,0.2256,0.2782,0.197,0.0925 0.3676,-0.1521,0.0591,0.0336,0.2227,0.2605,0.1985,0.0974 0.3824,-0.1591,0.0686,0.0293,0.2189,0.2415,0.1986,0.1005 0.3971,-0.1656,0.0786,0.0257,0.2143,0.2221,0.1971,0.1021 0.4118,-0.1714,0.0891,0.0229,0.209,0.2027,0.1941,0.1024 0.4265,-0.1765,0.0999,0.0209,0.2029,0.184,0.1896,0.1019 0.4412,-0.1811,0.111,0.02,0.1962,0.1664,0.1838,0.1009 0.4559,-0.1849,0.1222,0.0202,0.189,0.1505,0.1767,0.0998 0.4706,-0.1882,0.1335,0.0215,0.1812,0.1366,0.1685,0.099 0.4853,-0.1908,0.1447,0.024,0.1731,0.1251,0.1594,0.0987 0.5,-0.1928,0.1559,0.0278,0.1647,0.1163,0.1497,0.0993 0.5147,-0.1941,0.1668,0.0329,0.1562,0.1102,0.1396,0.101 0.5294,-0.1948,0.1774,0.0392,0.1475,0.1071,0.1292,0.104 0.5441,-0.1948,0.1877,0.0467,0.1388,0.1069,0.119,0.1082 0.5588,-0.1942,0.1974,0.0555,0.1303,0.1095,0.109,0.1138 0.5735,-0.193,0.2066,0.0655,0.1219,0.1149,0.0995,0.1206 0.5882,-0.1911,0.2151,0.0765,0.1138,0.1227,0.0906,0.1286 0.6029,-0.1886,0.2228,0.0886,0.1061,0.1328,0.0825,0.1374 0.6176,-0.1855,0.2297,0.1016,0.0988,0.1447,0.0752,0.1468 0.6324,-0.1817,0.2356,0.1154,0.092,0.1581,0.0688,0.1567 0.6471,-0.1773,0.2405,0.1299,0.0858,0.1724,0.0634,0.1666 0.6618,-0.1722,0.2443,0.1449,0.0802,0.1872,0.0587,0.1763 0.6765,-0.1665,0.2468,0.1602,0.0752,0.2019,0.0548,0.1854 0.6912,-0.1602,0.248,0.1758,0.0708,0.2159,0.0515,0.1939 0.7059,-0.1532,0.2479,0.1913,0.0671,0.2287,0.0485,0.2014 0.7206,-0.1456,0.2462,0.2066,0.0639,0.2397,0.0455,0.2078 0.7353,-0.1374,0.2429,0.2214,0.0613,0.2483,0.0423,0.2132 0.75,-0.1285,0.2379,0.2356,0.0592,0.254,0.0385,0.2176 0.7647,-0.1189,0.2312,0.2487,0.0574,0.2565,0.0337,0.2212 0.7794,-0.1088,0.2225,0.2606,0.056,0.2551,0.0276,0.2243 0.7941,-0.098,0.212,0.2709,0.0547,0.2498,0.0198,0.2272 0.8088,-0.0865,0.1993,0.2794,0.0533,0.2402,0.01,0.2304 0.8235,-0.0744,0.1845,0.2856,0.0517,0.2264,-0.0022,0.2345 0.8382,-0.0617,0.1675,0.2893,0.0497,0.2085,-0.0169,0.2399 0.8529,-0.0484,0.1482,0.2901,0.0469,0.1867,-0.0342,0.2471 0.8676,-0.0344,0.1264,0.2876,0.043,0.1615,-0.0542,0.2564 0.8824,-0.0197,0.1021,0.2814,0.0379,0.1337,-0.0765,0.2682 0.8971,-0.0045,0.0752,0.271,0.0309,0.1043,-0.1009,0.2821 0.9118,0.0115,0.0456,0.2562,0.0219,0.0746,-0.1265,0.2976 0.9265,0.028,0.0132,0.2363,0.0103,0.0462,-0.1525,0.3136 0.9412,0.0452,-0.0221,0.211,-0.0044,0.0212,-0.1774,0.3281 0.9559,0.063,-0.0604,0.1798,-0.0228,0.0018,-0.1994,0.3382 0.9706,0.0815,-0.1017,0.1421,-0.0453,-0.0091,-0.2164,0.3398 0.9853,0.1006,-0.1462,0.0975,-0.0726,-0.0083,-0.2255,0.3273
Gnuplotilla saa datan tulostettua viivoiksi:
plot for [col=2:8] 'data.csv' using 1:col with lines
Kun nyt kuitenkin haluat esittää, miten sikahienosti recomplex-luvuilla voi tehdä tämän pienimmän neliösumman approksimaation polynomeista, koodasin kyseisen menetelmän C++:lla ilman recomplex-lukuja tai valmiita matematiikkakirjastoja. En ole mikään matriisiekspertti ja tein tämän muiden hommien ohella parissa tunnissa, joten koodi ei varmaan ole täydellinen (enkä ole siinä viilannut koodaustyyliä todellakaan), mutta tuloksena on ihan oikeita polynomeja. Koodista pitää säätää datan mukaan toivottua polynomin astetta ja nollan toleranssia, ja myös polynomien määrä pitää tietää valmiiksi, kun en jaksanut tehdä parempaa CSV-datan käsittelyä. Koodi tulostaa Gnuplot-ohjelmalle sopivan plot-komennon, jolla voi piirtää kuvaajat.
Laitoin myös kuvan, josta voi todeta, että kohtuullisen samalta näyttää.
#include <fstream> #include <vector> #include <cassert> #include <cmath> #include <stdio.h> struct matriisi { const int w; int h; std::vector<double> data; matriisi(int w_, int h_): w(w_), h(h_), data(w * h) { } double& operator()(int x, int y) { return data[y * w + x]; } const double operator()(int x, int y) const { return data[y * w + x]; } }; matriisi matriisitulo(matriisi a, matriisi b) { assert(a.w == b.h || !"Tuloon tarvitaan yhteensopivat matriisit."); matriisi c(b.w, a.h); for (int x = 0; x < b.w; ++x) { for (int y = 0; y < a.h; ++y) { for (int k = 0; k < a.w; ++k) { c(x, y) += a(k, y) * b(x, k); } } } return c; } matriisi transpoosi(matriisi a) { matriisi b(a.h, a.w); for (int i = 0; i < a.w; ++i) { for (int j = 0; j < a.h; ++j) { b(j, i) = a(i, j); } } return b; } void plusrivi(matriisi& a, int kohde, int lahde, double kerroin) { for (int i = 0; i < a.w; ++i) { a(i, kohde) += kerroin * a(i, lahde); } } bool nolla(double d) { return std::abs(d) < 1.0e-10; } #define tulosta(m) tulosta_matriisi(m, # m) void tulosta_matriisi(matriisi const& a, const char* nimi = "") { printf("matriisi '%s' = %d x %d:\n", nimi, a.w, a.h); for (int y = 0; y < a.h; ++y) { for (int x = 0; x < a.w; ++x) { printf("%.4f, ", a(x, y)); } printf("\n"); } printf("\n"); } matriisi inversio(matriisi a) { assert(a.w == a.h || !"Inversioon tarvitaan neliömatriisi."); matriisi b(a.w, a.h); for (int i = 0; i < a.w; ++i) { b(i, i) = 1; } // Nollataan jokaisen rivin alkupuoli vähentämällä muita rivejä. for (int x = 0; x < a.h; ++x) { // Varmistetaan, että lävistäjällä on muu kuin nolla. if (nolla(a(x, x))) for (int y = x + 1; y < a.h; ++y) if (!nolla(a(x, y))) { double k = 1 / a(x, y); plusrivi(a, x, y, k); plusrivi(b, x, y, k); break; } if (nolla(a(x, x))) { throw std::runtime_error("Lävistäjälle jäi nolla, mitähän tälle tehdään?"); } // Skaalataan niin, että lävistäjällä on ykkönen. if (!nolla(a(x, x) - 1)) { double k = 1 / a(x, x) - 1; plusrivi(a, x, x, k); plusrivi(b, x, x, k); } // Vähennetään tätä riviä kaikista alemmista riveistä niin, // että alempien rivien y-sarakkeeseen jää pelkkää nollaa. for (int y = x + 1; y < a.h; ++y) if (a(x, y)) { double k = -a(x, y); plusrivi(a, y, x, k); plusrivi(b, y, x, k); } } // Nollataan jokaisen rivin loppupuoli vähentämällä muita rivejä. for (int x = a.h; x--;) { // Vähennetään tätä riviä kaikista ylemmista riveistä niin, // että ylempien rivien y-sarakkeeseen jää pelkkää nollaa. for (int y = x; y--;) if (a(x, y)) { double k = -a(x, y); plusrivi(a, y, x, k); plusrivi(b, y, x, k); } } return b; } int main() { auto polynomeja = 7; auto aste = 9; // x = datan ensimmäinen sarake ja sen potenssit. matriisi x(aste + 1, 0); // y = seuraavat sarakkeet. matriisi y(polynomeja, 0); // Luetaan syöte ja lasketaan x:n potenssit. // (Virheenkäsittely puuttuu, oletetaan oikean muotoista dataa.) std::ifstream ifs("data.csv"); if (!ifs) { throw std::runtime_error("data.csv puuttuu!"); } for (double luku; ifs >> luku;) { x.h += 1; x.data.push_back(1); x.data.push_back(luku); for (auto i = 2; i <= aste; ++i) { x.data.push_back(std::pow(luku, i)); } y.h += 1; for (int i = 0; i < polynomeja; ++i) { ifs.get(); ifs >> luku; y.data.push_back(luku); } } // Pienimmän neliösumman menetelmällä polynomin kertoimet ovat: // kertoimet = inversio(transpoosi(x) * x) * transpoosi(x) * y // B = inv(XT * X) * XT * y. matriisi xt = transpoosi(x); matriisi kertoimet = matriisitulo(matriisitulo(inversio(matriisitulo(xt, x)), xt), y); printf("gnuplot-komento:\n"); printf("plot [0:1] \\\n"); for (int i = 0; i < polynomeja; ++i) { printf(" %.3f + %.3f * x ", kertoimet(i, 0), kertoimet(i, 1)); for (int j = 2; j <= aste; ++j) { printf("+ %.3f * x**%d ", kertoimet(i, j), j); } puts((i+1 < polynomeja) ? " , \\" : ""); } #ifdef WIN32 printf("Enter lopettaa."), getchar(); #endif }
Mainittakoon vielä, että jollain valmiilla kirjastolla kuten Pythonin numpy tämän koodin saisi muutamaan riviin ja se toimisi vielä paremmin ja oudommissakin syötteissä.
Generoin vahingossa polynomit ilman kohinaa. Tai oikeastaan suodatin ne jo valmiiksi. Mutta ideana on, että input ja output korreloivat:
x=input[i]; y=poly.fx(x);
Käyrät ovat kvantisointuneet, jote kyllä niille periaatteessa pitää laskea pienimmän neliösumman approksimaatio.
Minun koodissani recomplex-oliopolynomi tulee lauseella:
polynomi P(input, output, 68/*pisteiden lkm*/, 8/*polynomin asteluku*/);
Aihe on jo aika vanha, joten et voi enää vastata siihen.