Sata lehtistipendiä jaettu
Tiede-lehti on jakanut reaaliaineissa menestyneille lukion oppilaille sata lehtistipendiä. Valitut saavat lehden vuosikerran. Stipendiaattien nimet löytyvät täältä.
|
|
KESKUSTELU
Tiede.fi-foorumin päävalikko. Keskustelua kaikille tieteestä kiinnostuneille. Edellyttää rekisteröitymistä.
Näytä vastaamattomat viestit | Näytä aktiiviset viestiketjut
| Kirjoittaja |
Viesti |
|
Läskiperse
|
 Lähetetty: Ti Heinä 10, 2012 3:28 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
|
Algebrallinen olio on ehkä vähiten väärä nimitys numeerisille olioille, jotka eivät oikein sovi kuntamäärittelyyn.
1) Jokaisella oliovariaatiolla on ensin alkeisoperaatiot plus-, miinus-, kerto- ja jakolasku sekä itseisarvo.
2) Oliot toteuttavat ehdot:
a*b = b*a
a/b=c <=> a=b*c
abs(a)*abs(b) = abs(a*b)
Nyt kaikki triviaalit ehdot pitävät myös automaattisesti paikkansa. Esim:
0+a = a+0
1*a = a
a+(-a) = 0
jne...
*** *** ***
Määritelmä pitää olla yksinkertainen ja efektiivinen. Matematiikassa kuitenkin tuhlataan rivejä itsestään selvien asioiden painamiseen.
Edellisestä oliovariaatioiden määrittelystä seuraa, että myös kaikki funktiot pitävät paikkansa, joihin päästään käsiksi sarjojen (alkeisoperaatioiden) kautta. Esim:
a = olio rnd()
b = sin(a)
c = arc sin(b)
Nyt a on varmasti yhtäsuuri kuin c.
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
 |
|
Läskiperse
|
Lähetetty: Ti Heinä 10, 2012 3:58 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
Oheinen koodi laskee 2D-reaalilukuina + 2D-kompleksioliolla. Koodi kääntyy takuulla jokaisessa C++ kääntäjässä ongelmitta. Lyhyellä main esimerkillä pääsee nopeasti kokeilemaan omia laskutuloksia. Cargo kirjoitti johonkin ketjuun negatiivisen logaritmin oton positiivisilla arvoilla. Jos sen paikan saisi uudelleen, voisin päivittää logaritmifunktiot. Koodi: #include <math.h>
#include <stdio.h>
class real2d
{
public:
double e[4];
real2d(void);
real2d(int);
real2d(int, int);
real2d(int, int, int, int);
real2d(double);
real2d(double, double);
friend int ifpos(real2d);
friend double abs(real2d);
friend void print(real2d);
friend real2d sqrt(real2d);
friend real2d reduse(real2d);
friend real2d operator-(real2d);
friend real2d operator+(real2d, real2d);
friend real2d operator-(real2d, real2d);
friend real2d operator*(real2d, real2d);
friend real2d operator/(real2d, real2d);
friend real2d operator*(double, real2d);
friend real2d operator/(real2d, double);
friend real2d operator*(real2d, double);
friend real2d operator/(double, real2d);
/* mgn = luvun suuruus */
friend double mgn(real2d);
friend real2d ln(real2d);
friend real2d exp(real2d);
friend real2d sin(real2d);
friend real2d cos(real2d);
friend real2d pow(real2d, int);
friend real2d pow(real2d, double);
friend real2d pow(real2d, real2d);
};
real2d::real2d(double a)
{
e[0]=a; e[1]=0;
e[2]=e[3]=0.00;
}
real2d::real2d(double a, double b)
{
e[0]=a; e[2]=b;
e[1]=e[3]=0.00;
}
real2d::real2d(int a)
{
*this=real2d((double)a);
}
real2d::real2d(int a, int b)
{
*this=real2d((double)a, (double)b);
}
real2d::real2d(int a, int b, int c, int d)
{
e[0]=a; e[1]=b;
e[2]=c; e[3]=d;
}
real2d::real2d(void)
{
*this=real2d((double)0);
}
real2d reduse(real2d a)
{
a.e[0]-=a.e[1];
a.e[2]-=a.e[3];
a.e[1]=a.e[3]=0;
return (real2d)a;
}
double abs(real2d a)
{
double x = a.e[0]-a.e[1];
double y = a.e[2]-a.e[3];
return (double)fabs(x-y);
}
double mgn(real2d a)
{
double x=a.e[0]-a.e[1];
double y=a.e[2]-a.e[3];
return fabs(x)+fabs(y);
}
int ifpos(real2d a)
{
double x=a.e[0]-a.e[1];
double y=a.e[2]-a.e[3];
if (x>=0.0)
{
return fabsl(y)<=x? 1: 0;
}
else
{
return 0;
}
}
real2d operator-(real2d a)
{
for (int i=0; i<4; i++)
a.e[i]=-a.e[i];
return reduse(a);
}
real2d operator+(real2d a, real2d b)
{
for (int i=0; i<4; i++)
a.e[i]+=b.e[i];
return reduse(a);
}
real2d operator-(real2d a, real2d b)
{
for (int i=0; i<4; i++)
a.e[i]-=b.e[i];
return reduse(a);
}
real2d operator*(real2d a, real2d b)
{
real2d t(0.0);
for (int i=0; i<4; i++)
for (int j=0; j<4; j++)
t.e[i^j]+=a.e[i]*b.e[j];
return reduse(t);
}
real2d operator/(real2d x, real2d y)
{
double X[4];
double K[4*4+4];
for (int n=0, p=0; n<4; n++)
{
for (int i=0; i<4; i++)
for (int j=0; j<4; j++)
{
if (n == (i^j))
{
K[p]=y.e[j]; ++p;
}
}
K[p]=x.e[n]; ++p;
}
// gsolve ratkaisee lineaarisen yhtälöryhmän
// Gaussin eliminointimenetelmällä...
int gsolve(double*, double*, int);
gsolve(X, K, 4);
double a=X[0]-X[1];
double b=X[2]-X[3];
return real2d(a, b);
}
void swapRows(double *a, double *b, int n)
{
for (int i=0; i<n; i++)
{
double tmp=a[i];
a[i]=b[i];
b[i]=tmp;
}
}
////////////////////////////////////////////////
// qsolve palautusarvo: //
// 1 <=> yhtälöryhmälle oli ratkaisu //
// 0 <=> yhtälöryhmälle ei ollut ratkaisua //
////////////////////////////////////////////////
int gsolve(double *X, double *f, int n)
{
double k, J;
double *max, *mux;
int i, j, x, M=(int)n+1;
for (i=0; i<n-1; i++)
{
max=f+i+i*M;
for (j=i+0x01; j < n; j++)
if (fabs(*(mux=f+i+j*M))>fabs(*max)) max=mux;
if (max!=(mux=f+i+i*M))
swapRows(max, mux, M-i);
if (fabs(J=f[i+i*M])<(double)1e-20)
return 0; /* eipä ole ratkaisuva */
for (j=i+1; j<n; j++)
{
if (f[i+j*M])
{
k=-f[i+j*M]/J;
for (x=i+1; x<M; x++)
f[x+j*M]=f[x+j*M]+k*f[x+i*M];
}
}
}
for (i=n-1; i>=0; i--)
{
k=(double)f[n+i*M];
for (j=n-1; j > i; j--)
k=k-(double)(X[j]*f[j+i*M]);
if (fabs(J=f[i+i*M])<(double)1e-20)
return 0; /* eipä ole ratkaisua */
else X[i]=k/J; /* onpa ratkaisu */
}
return 1;
}
real2d operator*(double k, real2d b)
{
real2d a(k);
return a*b;
}
real2d operator/(real2d a, double k)
{
real2d b(k);
return a/b;
}
real2d operator*(real2d a, double k)
{
real2d b(k);
return a*b;
}
real2d operator/(double k, real2d b)
{
real2d a(k);
return a/b;
}
real2d sqrt(real2d c)
{
if (ifpos(c))
{
real2d k(0.5);
real2d x(1.0);
for (int i=0; i<64; i++)
{
x=k*(x+c/x);
}
return x;
}
else
{
/* virhe, luku ei ole positiivinen */
return real2d(-(double)0xfedca201L);
}
}
void print(real2d a)
{
a = (real2d) reduse(a);
printf("(%0.21f, %0.21f)\n",
(double)a.e[0], (double)a.e[2]);
}
real2d pow(real2d x, int n)
{
if (n)
{
real2d t=x;
int i=n<0? -n: n;
for (--i; i; i--) t=t*x;
return n>0? t: real2d(1.0)/t;
}
else
{
return real2d(1.0);
}
}
real2d ln(real2d c)
{
int n=1;
double t;
real2d fx, dx;
real2d x1, x2;
real2d z(1.0);
real2d x(1.0), xx;
#define MAX 3.14159265e+36
if (!ifpos(c))
{
/* virhe, luku ei ole positiivinen */
return real2d(-(double)0xfedca201L);
}
while (abs(c)>1.001)
{
c=sqrt(c);
n+=n;
}
for (int k=2; k<160; k++)
for (int j=0; j<160; j++)
{
xx=x; t=1.0;
fx=z+x; dx=z;
for (int i=2; i<=k; i++)
{
dx=dx+i*xx/(t*=i);
xx=xx*x; fx=fx+xx/t;
if (mgn(fx)>MAX) break;
}
x1=x;
x=x-(fx-c)/dx;
x2=x;
if (k+1 < 160)
{
if (mgn(x1-x2)<1e-15) break;
}
}
return (double)n*x;
}
real2d exp(real2d x)
{
real2d xx;
double t=1.0;
real2d sum(1.0);
sum=sum+(xx=x);
for (int i=2; i<160; i++)
{
xx = xx*x;
sum=sum+xx/(t*=i);
if (mgn(xx)>MAX) break;
}
return sum;
}
real2d pow(real2d a, real2d x)
{
double t=1.0;
real2d sum(1.000);
real2d y=x*ln(a), xx=y;
sum=sum+xx;
for (int i=2; i<160; i++)
{
xx = xx*y;
sum=sum+xx/(t*=i);
if (mgn(xx)>MAX) break;
}
return sum;
}
real2d pow(real2d x, double ex)
{
return pow(x, real2d(ex));
}
real2d sin(real2d x)
{
double t=1.0;
real2d sum(0.0);
real2d y=x*x, xx=x;
for (int i=1, sg=0; i<160; i+=2, sg++)
{
if (sg&1) sum=sum-xx/t;
else /**/ sum=sum+xx/t;
xx=xx*y; t*=(i+1)*(i+2);
if (mgn(xx)>MAX) break;
}
return sum;
}
real2d cos(real2d x)
{
double t=2.0;
real2d sum(1.0);
real2d y=x*x, xx=y;
for (int i=2, sg=1; i<160; i+=2, sg++)
{
if (sg&1) sum=sum-xx/t;
else /**/ sum=sum+xx/t;
xx=xx*y; t*=(i+1)*(i+2);
if (mgn(xx)>MAX) break;
}
return sum;
}
class complex2d
{
public:
real2d r, i, j, k;
complex2d(void);
complex2d(int);
complex2d(double);
complex2d(real2d);
complex2d(double, double);
complex2d(real2d re, real2d im);
complex2d(real2d, real2d, real2d, real2d);
friend double abs(complex2d);
friend void print(complex2d);
friend complex2d sqrt(complex2d);
friend complex2d operator-(complex2d);
friend complex2d operator+(complex2d, complex2d);
friend complex2d operator-(complex2d, complex2d);
friend complex2d operator*(complex2d, complex2d);
friend complex2d operator/(complex2d, complex2d);
friend complex2d operator*(double, complex2d);
friend complex2d operator/(complex2d, double);
friend complex2d operator*(complex2d, double);
friend complex2d operator/(double, complex2d);
/* mgn = luvun suuruus */
friend double mgn(complex2d);
friend complex2d ln(complex2d);
friend complex2d exp(complex2d);
friend complex2d sin(complex2d);
friend complex2d cos(complex2d);
friend complex2d pow(complex2d, int);
friend complex2d pow(complex2d, double);
friend complex2d pow(complex2d, complex2d);
};
complex2d::complex2d(void)
{
r=i=j=k=real2d(0);
}
complex2d::complex2d(int a)
{
r=real2d(a, 0);
i=j=k=real2d(0);
}
complex2d::complex2d(double a)
{
r=real2d(a, 0.0);
i=j=k=real2d(0.0);
}
complex2d::complex2d(real2d a)
{
r=a; i=j=k=real2d(0);
}
complex2d::complex2d(double a, double b)
{
r=real2d(a, b);
i=j=k=real2d(0);
}
complex2d::complex2d(real2d a, real2d b)
{
r=a; i=b; j=k=real2d(0);
}
complex2d::complex2d(real2d a, real2d b, real2d c, real2d d)
{
r=a; i=b; j=c; k=d;
}
complex2d operator-(complex2d a)
{
return complex2d(-a.r, -a.i, -a.j, -a.k);
}
complex2d operator+(complex2d a, complex2d b)
{
return complex2d(a.r+b.r, a.i+b.i, a.j+b.j, a.k+b.k);
}
complex2d operator-(complex2d a, complex2d b)
{
return complex2d(a.r-b.r, a.i-b.i, a.j-b.j, a.k-b.k);
}
complex2d operator*(complex2d a, complex2d b)
{
complex2d c;
c.r = real2d(1,0,0,0)*a.r*b.r + real2d(0,1,0,0)*a.i*b.i
+ real2d(0,0,1,0)*a.j*b.j + real2d(0,0,0,1)*a.k*b.k;
c.i = real2d(1,0,0,0)*a.r*b.i + real2d(1,0,0,0)*a.i*b.r
+ real2d(1,0,0,0)*a.j*b.k + real2d(1,0,0,0)*a.k*b.j;
c.j = real2d(1,0,0,0)*a.r*b.j + real2d(0,0,0,1)*a.i*b.k
+ real2d(1,0,0,0)*a.j*b.r + real2d(0,0,0,1)*a.k*b.i;
c.k = real2d(1,0,0,0)*a.r*b.k + real2d(0,0,1,0)*a.i*b.j
+ real2d(0,0,1,0)*a.j*b.i + real2d(1,0,0,0)*a.k*b.r;
return c;
}
complex2d operator/(complex2d a, complex2d b)
{
complex2d L;
real2d z(0);
L=complex2d(b.r, b.i, -b.j, -b.k);
a = a*L;
b = b*L;
L=complex2d(b.r, -b.i, z, z);
a = a*L;
b = b*L;
a.r = a.r/b.r;
a.i = a.i/b.r;
a.j = a.j/b.r;
a.k = a.k/b.r;
return a;
}
complex2d operator*(double a, complex2d b)
{
return complex2d(a*b.r, a*b.i, a*b.j, a*b.k);
}
complex2d operator/(complex2d a, double b)
{
return complex2d(a.r/b, a.i/b, a.j/b, a.k/b);
}
complex2d operator*(complex2d a, double b)
{
return complex2d(b*a.r, b*a.i, b*a.j, b*a.k);
}
complex2d operator/(double a, complex2d b)
{
return complex2d(a)/b;
}
double mgn(complex2d a)
{
double sum=0;
sum+=mgn(a.r);
sum+=mgn(a.i);
sum+=mgn(a.j);
sum+=mgn(a.k);
return sum;
}
double abs(complex2d a)
{
complex2d L;
real2d z(0);
L=complex2d(a.r, a.i, -a.j, -a.k);
a = a*L;
L=complex2d(a.r, -a.i, z, z);
a = a*L;
a.r=pow(a.r, 0.25);
return abs(a.r);
}
complex2d sqrt(complex2d c)
{
if (mgn(c)<1e-20)
{
return complex2d(0);
}
else
{
complex2d x=c;
for (int i=0; i<160; i++)
x=x-(x*x-c)/(2.0*x);
return (complex2d)x;
}
}
complex2d sin(complex2d x)
{
double t=1.0;
complex2d sum(0.0);
complex2d y=x*x, xx=x;
for (int i=1, sg=0; i<160; i+=2, sg++)
{
if (sg&1) sum=sum-xx/t;
else /**/ sum=sum+xx/t;
xx=xx*y; t*=(i+1)*(i+2);
if (mgn(xx)>MAX) break;
}
return sum;
}
complex2d cos(complex2d x)
{
double t=2.0;
complex2d sum(1.0);
complex2d y=x*x, xx=y;
for (int i=2, sg=1; i<160; i+=2, sg++)
{
if (sg&1) sum=sum-xx/t;
else /**/ sum=sum+xx/t;
xx=xx*y; t*=(i+1)*(i+2);
if (mgn(xx)>MAX) break;
}
return sum;
}
complex2d pow(complex2d x, int n)
{
if (n)
{
complex2d t=x;
int i=n<0? -n: n;
for (--i; i; i--) t=t*x;
return n>0? t: complex2d(1.0)/t;
}
else
{
return complex2d(1.0);
}
}
complex2d ln(complex2d c)
{
int n=1;
double t;
complex2d fx, dx;
complex2d x1, x2;
complex2d z(1.0);
complex2d x(1.0), xx;
while (abs(c)>1.001)
{
c=sqrt(c);
n+=n;
}
for (int k=2; k<160; k++)
for (int j=0; j<160; j++)
{
xx=x; t=1.0;
fx=z+x; dx=z;
for (int i=2; i<=k; i++)
{
dx=dx+i*xx/(t*=i);
xx=xx*x; fx=fx+xx/t;
if (mgn(fx)>MAX) break;
}
x1=x;
x=x-(fx-c)/dx;
x2=x;
if (k+1 < 160)
{
if (mgn(x1-x2)<1e-15) break;
}
}
return (double)n*x;
}
complex2d exp(complex2d x)
{
complex2d xx;
double t=1.000;
complex2d sum(1);
sum=sum+(xx=x);
for (int i=2; i<160; i++)
{
xx = xx*x;
sum=sum+xx/(t*=i);
if (mgn(xx)>MAX) break;
}
return sum;
}
complex2d pow(complex2d a, complex2d x)
{
double t=(double)1;
complex2d sum(1.000);
complex2d ln(complex2d);
complex2d y=x*ln(a), xx=y;
sum=sum+xx;
for (int i=2; i<160; i++)
{
xx = xx*y;
sum=sum+xx/(t*=i);
if (mgn(xx)>MAX) break;
}
return sum;
}
complex2d pow(complex2d x, double n)
{
return pow(x, complex2d(n));
}
void print(complex2d a)
{
printf("[\n");
printf(" r="); print(a.r);
printf(" i="); print(a.i);
printf(" j="); print(a.j);
printf(" k="); print(a.k);
printf("]\n");
}
typedef unsigned long uint32;
uint32 rnd(uint32 max)
{
static uint32 R16A=1L;
static uint32 R16B=2L;
R16A -= 0x012357bfL;
R16B += 0xfedca201L;
R16A += (R16A>>16)^(R16B<<16);
R16B -= (R16A<<16)^(R16B>>16);
return (R16A^R16B)%max;
}
double drand(void)
{
uint32 jak=0xfffffffb;
double sg=rnd(2)? -1: 1;
double x=(double)rnd(jak);
return sg * x / (double)jak;
}
#pragma argsused
int main(int kpl, char* asia[])
{
for (int n=0; n<32767; n++)
{
real2d r(drand(), drand());
real2d i(drand(), drand());
real2d j(drand(), drand());
real2d k(drand(), drand());
complex2d c(r, i, j, k);
complex2d x=sqrt(c);
print(x);
print(x*x-c);
print(sin(c)*sin(c)+cos(c)*cos(c));
printf("================================\n\n");
}
return 0;
}
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Puuhikki
|
Lähetetty: Ti Heinä 10, 2012 5:10 pm |
|
Liittynyt: La Maalis 19, 2005 10:39 am
Viestit: 624
Paikkakunta: Kuopio
|
|
Ei se valitettavasti ihan noin mene. Siitä, että olet määritellyt magman ei seuraa, että olet määritellyt ykkösellisen magman. Ei ole myöskään triliaalia, että esim 0+a=0. Sinun pitää määritellä, että nolla on yhteenlaskun neutraalialkio.
En keksi, mitä hyötyä on määritellä joku random-algebrallinen struktuuri. Tietysti jos sillä pystyy lähestymään vaikkapa jotain otaksumaa, niin siitä vaan. Mutta tekstin perusteella en ole tästä kauhean vakuuttunut.
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: Ti Heinä 10, 2012 5:46 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
No ei se ihan randomilla ole määritelty, koska seula heittää kaikki epävalidit tapaukset saman tien pois. Sovelluksia olioille ei ole, paitsi polynomiapproksimaatiot. Nyt esim. real2d-oliolla voi laskea kahden käyrän pienimmän neliösumman polynomin virheen suhteessa vastinalkioihin, siis: X[0] = real2d(a0, b0), Y[0] = real2d(c0, 0) X[1] = real2d(a1, b1), Y[1] = real2d(c1, 0) X[2] = real2d(a2, b2), Y[2] = real2d(c2, 0) jne... => f(x) = a*x^0 + b*x^1 + c*x^2 ... Minä ole noita tutkinut koko elämäni. Ne oliot taipuvat ainakin korrelaatioiden laskemiseen. Senkin olen ajan kanssa huomannut, että aina pitää olla jokin 2^n akselinen viritys, että se pysyy mielekkäissä tuloksissa. Siis real2d:tä seuraa real4d, etc... Puuhikki kirjoitti: Ei ole myöskään triliaalia, että esim 0+a=0. Sinun pitää määritellä, että nolla on yhteenlaskun neutraalialkio. En oikein ymmärrä, mitä tarkoitat? Mutta triviaalisti määrittelystä seuraa automaattisesti, että: 0+a = asyy->seuraus = logiikka
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Puuhikki
|
Lähetetty: Ti Heinä 10, 2012 7:03 pm |
|
Liittynyt: La Maalis 19, 2005 10:39 am
Viestit: 624
Paikkakunta: Kuopio
|
|
Tarkoitan sitä, että nolla on määriteltävä jotenkin, ja yleensä tämä määritelmä on, että 0+a=a+0=a kaikilla a. Mutta koska 0+a=a seuraa jostakin, olisi kiva tietää miten määrittelet nollan.
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: Ti Heinä 10, 2012 7:50 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
|
Minun (seulan) ei tarvitse määritellä nollaa erikseen, koska kaikilla oliovariaatioilla on:
memset(*this, 0, sizeof(olio))
Olen myös pohtinut, että ehkä nämä oliot ovat sen vuoksi kummallisia, koska niiden kuvauskieli on enemmän jonkin ohjelmointikielen logiikkaa.
Sanoin aloituksessa, että "Algebrallinen olio on ehkä vähiten väärä nimitys numeerisille olioille."
Oliot itsessään sisältävät syy-seuraus määrityksiä sen jälkeen, kun on ensin asetettu oliolle ehdot:
a*b = b*a
a/b=c <=> a=b*c
abs(a)*abs(b) = abs(a*b)
Matemaattinen "tarpeellisuus" ilmenee vasta kontekstissa, eli kun taulukon syy-seurauksesta muodostetaan jonkin asteinen pienimmän neliösumman polynomi.
Yksi syy ehtojen efektisyydelle on tietokoneen teho.
Niin, ja sekin on oleellista, ettei noita olioita kukaan jaksa paperilla edes laskea.
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: Ke Heinä 11, 2012 2:47 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
Koodi: /******************************************************************************
Polynomiluokka ei sinällään ota kantaa sen käyttötarkoitukseen. *
Esim. kohinan poisto liukuvalla ikkunalla on toisinaan tarpeellista, etc. *
*
Kun koodin kääntää, sen käyttö on melko yksinkertaista. *
*
Ex. *
*
P = polynomi(nopeus, aika, näytepisteiden lukumäärä, polynomin asteluku) *
*
aika on taulukko: *
real2d aika[n]= *
{ *
real2d(t+0, 0), *
real2d(t+1, 0), *
real2d(t+2, 0), *
... *
}; *
*
nopeus on taulukko: *
real2d nopeus[n]= *
{ *
real2d(a0, b0), *
real2d(a1, b1), *
real2d(a2, b2), *
... *
}; *
*
näytepisteiden lukumäärä on yhtä kuin n *
*
polynomin asteluku voi olla joku 1, 2, 3, ..., eli ovat muotoa: *
*
1. aste: ax + b *
2. aste: ax^2 + bx + c *
3. aste: ax^3 + bx^2 + cx + d *
jne... *
*
P.fx(x) palauttaa funktion arvon arvolla x. y=P.fx(x) *
*
P.derivaatta(x) palauttaa derivaatan kohdassa x. *
*
P.integraali(a, b) palauttaa integraalin väliltä [a, b]. *
*
Kuta enemmän näytepisteitä tarkasteltavalta väliltä on, *
sen paremmin laskettu polynomi pistejoukkoa kuvaa. *
*
Kysymyksiä jos tulee, niihin on helpompi vastata, *
kun on kysymys, johon vastata. *
*
******************************************************************************/
Nyt joku voisi pyöräyttää jollakin ex. MatCad:lla kaksi nopeutta suhteessa lineaariseen aikaan. Kirjoitan sitten viimeisen funktion niillä arvoilla ja lasken polynomin. Sen jälkeen tulostan ko. aika-arvoilla nopeudet. Oheisesta koodista on poistettu kaikki ylimääräinen, joten kynnys viimeisen main-funktion kirjoittamiselle itse ei pitäisi olla suuri. Koodi: #include <math.h>
#include <stdio.h>
class real2d
{
public:
double e[4];
real2d(int);
real2d(int, int);
real2d(void);
real2d(double);
real2d(double, double);
friend int ifpos(real2d);
friend double abs(real2d);
friend void print(real2d);
friend real2d reduse(real2d);
friend real2d operator-(real2d);
friend real2d operator+(real2d, real2d);
friend real2d operator-(real2d, real2d);
friend real2d operator*(real2d, real2d);
friend real2d operator/(real2d, real2d);
friend real2d operator*(double, real2d);
friend real2d operator/(real2d, double);
friend real2d operator*(real2d, double);
friend real2d operator/(double, real2d);
/* mgn = luvun suuruus */
friend double mgn(real2d);
friend real2d pow(real2d, int);
};
real2d::real2d(double a)
{
e[0]=a; e[1]=0;
e[2]=e[3]=0.00;
}
real2d::real2d(double a, double b)
{
e[0]=a; e[2]=b;
e[1]=e[3]=0.00;
}
real2d::real2d(int a)
{
*this=real2d((double)a);
}
real2d::real2d(int a, int b)
{
*this=real2d((double)a, (double)b);
}
real2d::real2d(void)
{
*this=real2d((double)0);
}
real2d reduse(real2d a)
{
a.e[0]-=a.e[1];
a.e[2]-=a.e[3];
a.e[1]=a.e[3]=0;
return (real2d)a;
}
double abs(real2d a)
{
double x = a.e[0]-a.e[1];
double y = a.e[2]-a.e[3];
return (double)fabs(x-y);
}
double mgn(real2d a)
{
double x=a.e[0]-a.e[1];
double y=a.e[2]-a.e[3];
return fabs(x)+fabs(y);
}
int ifpos(real2d a)
{
double x=a.e[0]-a.e[1];
double y=a.e[2]-a.e[3];
if (x>=0.0)
{
return fabs(y)<=x? 1: 0;
}
else
{
return 0;
}
}
real2d operator-(real2d a)
{
for (int i=0; i<4; i++)
a.e[i]=-a.e[i];
return reduse(a);
}
real2d operator+(real2d a, real2d b)
{
for (int i=0; i<4; i++)
a.e[i]+=b.e[i];
return reduse(a);
}
real2d operator-(real2d a, real2d b)
{
for (int i=0; i<4; i++)
a.e[i]-=b.e[i];
return reduse(a);
}
real2d operator*(real2d a, real2d b)
{
real2d t(0.0);
for (int i=0; i<4; i++)
for (int j=0; j<4; j++)
t.e[i^j]+=a.e[i]*b.e[j];
return reduse(t);
}
real2d operator/(real2d x, real2d y)
{
double X[4];
double K[4*4+4];
for (int n=0, p=0; n<4; n++)
{
for (int i=0; i<4; i++)
for (int j=0; j<4; j++)
{
if (n == (i^j))
{
K[p]=y.e[j]; ++p;
}
}
K[p]=x.e[n]; ++p;
}
// gsolve ratkaisee lineaarisen yhtälöryhmän
// Gaussin eliminointimenetelmällä...
int gsolve(double*, double*, int);
gsolve(X, K, 4);
double a=X[0]-X[1];
double b=X[2]-X[3];
return real2d(a, b);
}
void swapRows(double *a, double *b, int n)
{
for (int i=0; i<n; i++)
{
double tmp=a[i];
a[i]=b[i];
b[i]=tmp;
}
}
////////////////////////////////////////////////
// qsolve palautusarvo: //
// 1 <=> yhtälöryhmälle oli ratkaisu //
// 0 <=> yhtälöryhmälle ei ollut ratkaisua //
////////////////////////////////////////////////
int gsolve(double *X, double *f, int n)
{
double k, J;
double *max, *mux;
int i, j, x, M=(int)n+1;
for (i=0; i<n-1; i++)
{
max=f+i+i*M;
for (j=i+0x01; j < n; j++)
if (fabs(*(mux=f+i+j*M))>fabs(*max)) max=mux;
if (max!=(mux=f+i+i*M))
swapRows(max, mux, M-i);
if (fabs(J=f[i+i*M])<(double)1e-20)
return 0; /* eipä ole ratkaisuva */
for (j=i+1; j<n; j++)
{
if (f[i+j*M])
{
k=-f[i+j*M]/J;
for (x=i+1; x<M; x++)
f[x+j*M]=f[x+j*M]+k*f[x+i*M];
}
}
}
for (i=n-1; i>=0; i--)
{
k=(double)f[n+i*M];
for (j=n-1; j > i; j--)
k=k-(double)(X[j]*f[j+i*M]);
if (fabs(J=f[i+i*M])<(double)1e-20)
return 0; /* eipä ole ratkaisua */
else X[i]=k/J; /* onpa ratkaisu */
}
return 1;
}
real2d operator*(double k, real2d b)
{
real2d a(k);
return a*b;
}
real2d operator/(real2d a, double k)
{
real2d b(k);
return a/b;
}
real2d operator*(real2d a, double k)
{
real2d b(k);
return a*b;
}
real2d operator/(double k, real2d b)
{
real2d a(k);
return a/b;
}
void print(real2d a)
{
a = (real2d) reduse(a);
printf("(%0.21f, %0.21f)\n",
(double)a.e[0], (double)a.e[2]);
}
real2d pow(real2d x, int n)
{
if (n)
{
real2d t=x;
int i=n<0? -n: n;
for (--i; i; i--) t=t*x;
return n>0? t: real2d(1.0)/t;
}
else
{
return real2d(1.0);
}
}
/* ** *** *** *** *** *** *** *** *** *** *** *** *** *** *** *** *** *** ** */
/* ** *** *** *** *** *** *** *** *** *** *** *** *** *** *** *** *** *** ** */
/* ** *** *** *** *** *** *** *** *** *** *** *** *** *** *** *** *** *** ** */
/******************************************************************************
Polynomiluokka ei sinällään ota kantaa sen käyttötarkoitukseen. *
Esim. kohinan poisto liukuvalla ikkunalla on toisinaan tarpeellista, etc. *
*
Kun koodin kääntää, sen käyttö on melko yksinkertaista. *
*
Ex. *
*
P = polynomi(nopeus, aika, näytepisteiden lukumäärä, polynomin asteluku) *
*
aika on taulukko: *
real2d aika[n]= *
{ *
real2d(t+0, 0), *
real2d(t+1, 0), *
real2d(t+2, 0), *
... *
}; *
*
nopeus on taulukko: *
real2d nopeus[n]= *
{ *
real2d(a0, b0), *
real2d(a1, b1), *
real2d(a2, b2), *
... *
}; *
*
näytepisteiden lukumäärä on yhtä kuin n *
*
polynomin asteluku voi olla joku 1, 2, 3, ..., eli ovat muotoa: *
*
1. aste: ax + b *
2. aste: ax^2 + bx + c *
3. aste: ax^3 + bx^2 + cx + d *
jne... *
*
P.fx(x) palauttaa funktion arvon arvolla x. y=P.fx(x) *
*
P.derivaatta(x) palauttaa derivaatan kohdassa x. *
*
P.integraali(a, b) palauttaa integraalin väliltä [a, b]. *
*
Kuta enemmän näytepisteitä tarkasteltavalta väliltä on, *
sen paremmin laskettu polynomi pistejoukkoa kuvaa. *
*
Kysymyksiä jos tulee, niihin on helpompi vastata, *
kun on kysymys, johon vastata. *
*
******************************************************************************/
class polynomi
{
public: /* olion rajapinta */
~polynomi(void);
polynomi(real2d*, real2d*, int, int);
real2d fx(real2d x);
real2d derivaatta(real2d);
real2d integraali(real2d, real2d);
private: /* olion sisäiset parametrit ja funktiot */
int p_aste;
real2d *K;
double gabs(real2d);
void mem00set(void*, int);
real2d inteGraali(real2d);
real2d getSumX(real2d*, int);
real2d getSumXY(real2d*, real2d*, int);
real2d fxn(real2d*, int, real2d);
int ratkaise(real2d*, real2d*, int);
void swap_rows(real2d*, real2d*, int);
void mk_next_pow(real2d*, real2d*, int);
void mk_poly(real2d*, real2d*, int, int);
void copyMatrixLeft(real2d*, real2d, int, int);
void copyMatrixDown(real2d*, real2d, int, int);
};
polynomi::polynomi(real2d *x, real2d *y, int n, int a)
{
p_aste=a;
K=new real2d[a+1];
mk_poly(x, y, n, p_aste);
}
polynomi::~polynomi(void)
{
delete[] K;
}
real2d polynomi::fx(real2d x)
{
return fxn(K, p_aste, x);
}
real2d polynomi::derivaatta(real2d x)
{
real2d xxx(1.0);
real2d sum=real2d(0.0);
for (int i=1; i<=p_aste; i++)
{
sum=sum+K[i]*real2d(i)*xxx;
xxx=xxx*x;
}
return sum;
}
real2d polynomi::inteGraali(real2d x)
{
real2d xxx=x*x;
real2d sum=real2d(0);
for (int i=1; i<=p_aste; i++)
{
real2d z=real2d(1)/real2d(i+1);
sum=sum+z*K[i]*xxx;
xxx=xxx*x;
}
return sum;
}
real2d polynomi::integraali(real2d a, real2d b)
{
return inteGraali(b)-inteGraali(a);
}
double polynomi::gabs(real2d x)
{
return mgn(x);
}
void polynomi::mem00set(void *x, int kpl)
{
char *c=(char*)x;
for (int i=0; i<kpl; i++)
c[i]=(char)0;
}
void polynomi::swap_rows(real2d *a, real2d *b, int n)
{
for (int i=0; i<n; i++)
{
real2d tmp=a[i];
a[i]=b[i];
b[i]=tmp;
}
}
int polynomi::ratkaise(real2d *X, real2d *f, int n)
{
real2d k, J;
real2d *max, *mux;
int i, j, x, M=(int)n+1;
for (i=0; i<n-1; i++)
{
max=f+i+i*M;
for (j=i+1; j<n; j++)
if (gabs(*(mux=f+i+j*M))>gabs(*max)) max=mux;
if (max!=(mux=f+i+i*M))
swap_rows(max, mux, M-i);
if (gabs(J=f[i+i*M])<(double)1e-15)
return 0; /* eipä ole ratkaisua */
for (j=i+1; j<n; j++)
{
if (gabs(f[i+j*M]))
{
k=-f[i+j*M]/J;
for (x=i+1; x<M; x++)
f[x+j*M]=f[x+j*M]+k*f[x+i*M];
}
}
}
for (i=n-1; i>=0; i--)
{
k=(real2d)f[n+i*M];
for (j=n-1; j > i; j--)
k=k-(real2d)(X[j]*f[j+i*M]);
if (gabs(J=f[i+i*M])<(double)1e-15)
return 0; /* eipä ole ratkaisua */
else X[i]=k/J; /* onpa ratkaisu */
}
return 1;
}
real2d polynomi::fxn(real2d *k, int n, real2d x)
{
real2d sum=k[0];
real2d y=real2d(1);
for (int i=1; i<=n; i++)
{
y=y*x;
sum=sum+k[i]*y;
}
return sum;
}
real2d polynomi::getSumX(real2d *x, int n)
{
real2d sum=real2d(0);
for (int i=0; i<n; i++)
sum=sum+x[i]; return sum;
}
real2d polynomi::getSumXY(real2d *x, real2d *y, int n)
{
real2d sum=real2d(0);
for (int i=0; i<n; i++)
sum=sum+x[i]*y[i]; return sum;
}
void polynomi::mk_next_pow(real2d *x, real2d *y, int n)
{
for (int i=0; i<n; i++) x[i]=x[i]*y[i];
}
void polynomi::copyMatrixLeft(real2d *L, real2d value, int i, int N)
{
for (--i; i>=0; i--) {L+=N; *L=value;}
}
void polynomi::copyMatrixDown(real2d *L, real2d value, int i, int N)
{
for (i+=1; i<N; i++) {L+=N; *L=value;}
}
void polynomi::mk_poly(real2d *X, real2d *Y, int kpl, int aste)
{
int N=aste+1, i;
int M=aste+2, j;
int upa = N*N+N;
real2d *x=new real2d[kpl];
real2d *L=new real2d[upa];
mem00set(x, sizeof(real2d)*kpl);
mem00set(L, sizeof(real2d)*upa);
for (i=0; i<kpl; i++) x[i]=(real2d)1;
L[N]=getSumXY(x, Y, (int)kpl);
L[0]=getSumX(x, kpl);
for (i=1; i<N; i++)
{
mk_next_pow(x, X, kpl);
L[i]=getSumX(x, (int)kpl);
L[M*i+N]=getSumXY(x, Y, kpl);
copyMatrixLeft(L+i, L[i], i, N);
}
for (i=1; i<N; i++)
{
j=(int)(i*M+N-0x01);
mk_next_pow(x, X, kpl);
L[j]=getSumX(x, (int)kpl);
copyMatrixDown(L+j, L[j], i, N);
}
if (!ratkaise(K, L, N))
{
/* eipä ollut ratkaisua */
mem00set(K, sizeof(real2d)*N);
}
delete[] L;
delete[] x;
}
#pragma argsused
int main(int kpl, char* asia[])
{
return 0;
}
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: Ke Heinä 11, 2012 9:19 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
Eräs ilmiö liittyi aikaan siten, että ajanhetkellä 0 mitattiin nopeus(14, 8). Ajanhetkellä 1 mitattiin nopeus(-3, 10). Todettiin, että ilmiön syy-seuraus noudatti likimain arvoja: aika0=(0, 0), nopeus0=(14, 8) aika1=(1, 0), nopeus1=(-3, 10) aika2=(2, 0), nopeus2=(-12, 14) aika3=(3, 0), nopeus3=(-20, 17) aika4=(4, 0), nopeus4=(-25, 19) aika5=(5, 0), nopeus5=(-28, 22) aika6=(6, 0), nopeus6=(-27, 24) aika7=(7, 0), nopeus7=(-25, 27) aika8=(8, 0), nopeus8=(-23, 27) aika9=(9, 0), nopeus9=(-18, 27) aika10=(10, 0), nopeus10=(-13, 28) aika11=(11, 0), nopeus11=(-4, 28) aika12=(12, 0), nopeus12=(1, 28) aika13=(13, 0), nopeus13=(9, 27) aika14=(14, 0), nopeus14=(15, 26) aika15=(15, 0), nopeus15=(22, 24) aika16=(16, 0), nopeus16=(26, 21) aika17=(17, 0), nopeus17=(30, 21) aika18=(18, 0), nopeus18=(33, 17) aika19=(19, 0), nopeus19=(34, 13) aika20=(20, 0), nopeus20=(32, 11) aika21=(21, 0), nopeus21=(28, 7) aika22=(22, 0), nopeus22=(22, 2) aika23=(23, 0), nopeus23=(12, -2) aika24=(24, 0), nopeus24=(-1, -8) (Kohina on lisätty randomilla.) Graafisen tarkastelun perusteella ilmiölle päätettiin sovittaa kuutioparaabeli. Koodi: #pragma argsused
int main(int kpl, char* asia[])
{
real2d aika[25];
real2d nopeus[25];
aika[ 0]=real2d( 0, 0), nopeus[ 0]=real2d( 14, 8);
aika[ 1]=real2d( 1, 0), nopeus[ 1]=real2d( -3, 10);
aika[ 2]=real2d( 2, 0), nopeus[ 2]=real2d(-12, 14);
aika[ 3]=real2d( 3, 0), nopeus[ 3]=real2d(-20, 17);
aika[ 4]=real2d( 4, 0), nopeus[ 4]=real2d(-25, 19);
aika[ 5]=real2d( 5, 0), nopeus[ 5]=real2d(-28, 22);
aika[ 6]=real2d( 6, 0), nopeus[ 6]=real2d(-27, 24);
aika[ 7]=real2d( 7, 0), nopeus[ 7]=real2d(-25, 27);
aika[ 8]=real2d( 8, 0), nopeus[ 8]=real2d(-23, 27);
aika[ 9]=real2d( 9, 0), nopeus[ 9]=real2d(-18, 27);
aika[10]=real2d(10, 0), nopeus[10]=real2d(-13, 28);
aika[11]=real2d(11, 0), nopeus[11]=real2d( -4, 28);
aika[12]=real2d(12, 0), nopeus[12]=real2d( 1, 28);
aika[13]=real2d(13, 0), nopeus[13]=real2d( 9, 27);
aika[14]=real2d(14, 0), nopeus[14]=real2d( 15, 26);
aika[15]=real2d(15, 0), nopeus[15]=real2d( 22, 24);
aika[16]=real2d(16, 0), nopeus[16]=real2d( 26, 21);
aika[17]=real2d(17, 0), nopeus[17]=real2d( 30, 21);
aika[18]=real2d(18, 0), nopeus[18]=real2d( 33, 17);
aika[19]=real2d(19, 0), nopeus[19]=real2d( 34, 13);
aika[20]=real2d(20, 0), nopeus[20]=real2d( 32, 11);
aika[21]=real2d(21, 0), nopeus[21]=real2d( 28, 7);
aika[22]=real2d(22, 0), nopeus[22]=real2d( 22, 2);
aika[23]=real2d(23, 0), nopeus[23]=real2d( 12, -2);
aika[24]=real2d(24, 0), nopeus[24]=real2d( -1, -8);
polynomi P(aika, nopeus, 25, 3);
for (int i=0; i<25; i++)
{
real2d noPeus = P.fx(aika[i]);
print(noPeus);
}
return 0;
} Kun mittausvastinarvot sijoitetaan pienimmän neliösumman kaavaan, kuutioparaabelin antamat nopeudet (P.fx(aika_n)) ovat: Koodi: ( 13.131, 7.000)
( -1.394, 10.758)
(-12.409, 14.152)
(-20.227, 17.181)
(-25.163, 19.842)
(-27.528, 22.134)
(-27.638, 24.053)
(-25.805, 25.598)
(-22.342, 26.766)
(-17.563, 27.554)
(-11.782, 27.961)
( -5.311, 27.984)
( 1.536, 27.621)
( 8.445, 26.870)
( 15.103, 25.728)
( 21.197, 24.192)
( 26.413, 22.261)
( 30.438, 19.932)
( 32.959, 17.203)
( 33.661, 14.072)
( 32.232, 10.535)
( 28.359, 6.592)
( 21.727, 2.239)
( 12.024, -2.525)
( -1.064, -7.704)
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: To Heinä 12, 2012 10:14 am |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
|
Jos edellisen main:in alkuperäiset seuraukset piirtää kuutioparaabelin antamien seurauksien keralla, kuvasta näkee, kuinka kohina on eliminoitunut.
Jos mittausvastinarvojen syy ja seuraus kulkevat toisiaan vasten, polynomista voi tulla myös kompleksinen.
Jotta esitys olisi yleisempi, mittausvastinparien nimet olkoot syy[n] ja seuraus[n]. Nyt päästään epäsymmetrisiin syy-seurauksiin.
Esim. real4d:llä voi olla korrelaatio:
syy0=(a0, b0, 0, 0), seuraus0=(c0, d0, e0, f0)
syy1=(a1, b1, 0, 0), seuraus1=(c1, d1, e1, f1)
syy2=(a2, b2, 0, 0), seuraus2=(c2, d2, e2, f2)
...
Tai real2d:llä voi olla korrelaatio:
syy0=(a0, b0), seuraus0=(c0, d0)
syy1=(a1, b1), seuraus1=(c1, d1)
syy2=(a2, b2), seuraus2=(c2, d2)
...
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: To Heinä 12, 2012 1:44 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
|
real2d:ssä on neljä alkiota (double e[4]). Syy siihen on, että kertolasku ja jakolasku ovat mahdollista suorittaa vain redusoimattomana. Parittomat alkiot ovat parillisien alkioiden apu-arvoja.
Apu-arvoista pääsi eroon, kun vain kerto- ja jakolasku suoritetaan redusoimattomana varaamalla niihin funktioiden sisällä kaksinkertaiset taulukot.
real4d:stä tuli puolet lyhyempi. Ja varsinainen pointti matematiikan professoreille on, että eikö ole koskaan tullut mieleen suorittaa alkeisoperaatioita redusoimattomana?
Aikamoinen suo on jäänyt korrelaatioiden suhteen kyntämättä. Ja siksi heristän nyrkkiäni niin usein matematiikalle.
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: To Heinä 12, 2012 3:11 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
Oheisen real4d-koodin viimeinen main-funktio testaa alkeisoperaatioita aloitusviestin mukaan. Kerto- ja jakolasku suoritetaan ensin redusoimattomana, ja funktioiden lopussa suoritetaan redusointi. Koodi: #include <math.h>
#include <stdio.h>
class real4d
{
public:
double f[4];
real4d(int);
real4d(int, int);
real4d(int, int, int, int);
real4d(void);
real4d(double);
real4d(double, double);
real4d(double, double, double, double);
friend double abs(real4d);
friend void print(real4d);
friend real4d operator-(real4d);
friend real4d operator+(real4d, real4d);
friend real4d operator-(real4d, real4d);
friend real4d operator*(real4d, real4d);
friend real4d operator*(double, real4d);
friend real4d operator*(real4d, double);
friend real4d operator/(real4d, real4d);
friend real4d operator/(double, real4d);
friend real4d operator/(real4d, double);
/* mgn=luvun suuruus */
friend double mgn(real4d c);
friend void mem00set(void*, int);
};
void mem00set(void *ptr, int kpl)
{
char *c=(char*)ptr;
for (int i=0; i<kpl; i++)
c[i]=(char)0;
}
real4d::real4d(void)
{
mem00set(this, sizeof(real4d));
}
real4d::real4d(int a)
{
mem00set(this, sizeof(real4d));
f[0]=a;
}
real4d::real4d(int a, int b)
{
mem00set(this, sizeof(real4d));
f[0]=a; f[1]=b;
}
real4d::real4d(int a, int b, int c, int d)
{
mem00set(this, sizeof(real4d));
f[0]=a; f[1]=b; f[2]=c; f[3]=d;
}
real4d::real4d(double a)
{
mem00set(this, sizeof(real4d));
f[0]=a;
}
real4d::real4d(double a, double b)
{
mem00set(this, sizeof(real4d));
f[0]=a; f[1]=b;
}
real4d::real4d(double a, double b, double c, double d)
{
mem00set(this, sizeof(real4d));
f[0]=a; f[1]=b; f[2]=c; f[3]=d;
}
double abs(real4d c)
{
double *f=c.f;
return fabs(f[0]-f[1]+f[2]-f[3]);
}
double mgn(real4d c)
{
double sum=0;
for (int i=0; i<4; i++)
sum+=fabs(c.f[i]);
return sum;
}
real4d operator-(real4d a)
{
for (int i=0; i<4; i++)
a.f[i]=-a.f[i];
return a;
}
real4d operator+(real4d a, real4d b)
{
for (int i=0; i<4; i++)
a.f[i]+=b.f[i];
return a;
}
real4d operator-(real4d a, real4d b)
{
for (int i=0; i<4; i++)
a.f[i]-=b.f[i];
return a;
}
real4d operator*(real4d a, real4d b)
{
double x[8];
double y[8];
double t[8];
mem00set(x, sizeof(double)*8);
mem00set(y, sizeof(double)*8);
mem00set(t, sizeof(double)*8);
for (int i=0; i<4; i++)
{
x[i*2]=a.f[i];
y[i*2]=b.f[i];
}
for (int i=0; i<8; i++)
for (int j=0; j<8; j++)
t[i^j]+=x[i]*y[j];
for (int i=0; i<4; i++)
{
t[i]=t[i*2]-t[i*2+1];
}
return real4d(t[0], t[1], t[2], t[3]);
}
real4d operator*(double k, real4d b)
{
real4d a(k);
return a*b;
}
real4d operator*(real4d a, double k)
{
real4d b(k);
return a*b;
}
real4d operator/(double k, real4d b)
{
real4d a(k);
return a/b;
}
real4d operator/(real4d a, double k)
{
real4d b(k);
return a/b;
}
real4d operator/(real4d u, real4d v)
{
double x[8];
double y[8];
mem00set(x, sizeof(double)*8);
mem00set(y, sizeof(double)*8);
for (int i=0; i<4; i++)
{
x[i*2]=u.f[i];
y[i*2]=v.f[i];
}
double X[8];
double K[8*8+8];
for (int n=0, p=0; n<8; n++)
{
for (int i=0; i<8; i++)
for (int j=0; j<8; j++)
{
if (n == (i^j))
{
K[p]=y[j]; ++p;
}
}
K[p]=x[n]; ++p;
}
// gsolve ratkaisee lineaarisen yhtälöryhmän
// Gaussin eliminointimenetelmällä...
int gsolve(double*, double*, int);
gsolve(X, K, 8);
double a=X[0]-X[1];
double b=X[2]-X[3];
double c=X[4]-X[5];
double d=X[6]-X[7];
return real4d(a, b, c, d);
}
////////////////////////////////////////////////
// qsolve palautusarvo: //
// 1 <=> yhtälöryhmälle oli ratkaisu //
// 0 <=> yhtälöryhmälle ei ollut ratkaisua //
////////////////////////////////////////////////
int gsolve(double *X, double *f, int n)
{
double k, J;
double *max, *mux;
int i, j, x, M=(int)n+1;
void swapRows(double*, double*, int);
for (i=0; i<n-1; i++)
{
max=f+i+i*M;
for (j=i+0x01; j < n; j++)
if (fabs(*(mux=f+i+j*M))>fabs(*max)) max=mux;
if (max!=(mux=f+i+i*M))
swapRows(max, mux, M-i);
if (fabs(J=f[i+i*M])<(double)1e-15)
return 0; /* eipä ole ratkaisuva */
for (j=i+1; j<n; j++)
{
if (f[i+j*M])
{
k=-f[i+j*M]/J;
for (x=i+1; x<M; x++)
f[x+j*M]=f[x+j*M]+k*f[x+i*M];
}
}
}
for (i=n-1; i>=0; i--)
{
k=(double)f[n+i*M];
for (j=n-1; j > i; j--)
k=k-(double)(X[j]*f[j+i*M]);
if (fabs(J=f[i+i*M])<(double)1e-15)
return 0; /* eipä ole ratkaisua */
else X[i]=k/J; /* onpa ratkaisu */
}
return 1;
}
void swapRows(double *a, double *b, int n)
{
for (int i=0; i<n; i++)
{
double tmp=a[i];
a[i] = b[i];
b[i]=tmp;
}
}
void print(real4d a)
{
printf("(%0.10f, %0.10f, %0.10f, %0.10f)\n",
(double)a.f[0], (double)a.f[1], (double)a.f[2], (double)a.f[3]);
}
/* ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** */
/* ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** */
/* ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** */
// tuotetaan satunnaisluku
unsigned long rnd(unsigned long max)
{
static unsigned long R16A=1L;
static unsigned long R16B=2L;
R16A -= 0x012357bfL;
R16B += 0xfedca201L;
R16A += (R16A>>16)^(R16B<<16);
R16B -= (R16A<<16)^(R16B>>16);
return (R16A^R16B)%max;
}
// tuotetaan satunnaisluku
double drand(void)
{
unsigned long jak=0xfffffffb;
double sg=rnd(0x02)? -1: 1;
double x=(double)rnd(jak);
return sg * x/(double)jak;
}
#pragma argsused
int main(int kpl, char* asia[])
{
for (int i=0; i<16; i++)
{
real4d a(drand(), drand(), drand(), drand());
real4d b(drand(), drand(), drand(), drand());
printf("ehto 1: a*b = b*a \n");
print(a*b);
print(b*a);
printf("\n");
printf("ehto 2: a/b=c <=> a=b*c \n");
real4d c=a/b;
print(a);
print(b*c);
printf("\n");
printf("ehto 3: abs(a)*abs(b) = abs(a*b) \n");
printf("%25.20f\n", abs(a)*abs(b));
printf("%25.20f\n", abs(a*b));
printf("\n");
}
return 0;
}
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: To Heinä 12, 2012 6:25 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
huoh...pitää varmasti seuraavaksi tehdä yleis-versio, joka tukee kaikkia ulottuvuuksia, siis jokin: Koodi: ...
#define DIM (1<<n);
double f[DIM];
... jossa n voi olla 0, 1, 2, 4, ... nythän matematiikka on suuressa viisaudessaan päätynyt siihen, että tässä meillä nyt on yksiulotteinen reaaliluku: Koodi: ...
#define DIM (1<<n);
double f[DIM];
... jossa n voi olla vain ja ainoastaan 0! väännän mukaan polynomiapproksimaatiot reaalisissa ja kompleksisissa sfääreissä. palaan astialle...
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: To Heinä 12, 2012 8:32 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
tai saman asian voi sanoa paljon selkeämminkin (ilman bittisiirtoja), jossa korostuu sarja 2^n. seuraavaksi on tehdä yleis-versio, joka tukee kaikkia ulottuvuuksia, mutta siis: Koodi: ...
#define DIM 4
double f[DIM];
... jossa DIM voi olla 1, 2, 4, 8, ... nythän matematiikka on suuressa viisaudessaan päätynyt siihen, että tässä meillä nyt on yksiulotteinen reaaliluku: Koodi: ...
#define DIM 1
double f[DIM];
... jossa DIM voi olla vain ja ainoastaan 1!
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: To Heinä 12, 2012 11:09 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
ohessa on yleis areal aloitusviestissä mainittujen ehtojen testaamiseen.
areal:ssa DIM voi olla 1, 2, 4, 8, ...Koodi: #include <stdio.h>
class areal
{
public:
#define DIM 8
double f[DIM];
areal(void);
areal(double*);
areal(char*, ...);
friend void print(areal);
friend double abs(areal);
friend double abs(double);
friend areal operator-(areal);
friend areal operator+(areal, areal);
friend areal operator-(areal, areal);
friend areal operator*(areal, areal);
friend areal operator*(double, areal);
friend areal operator*(areal, double);
friend areal operator/(areal, areal);
friend areal operator/(double, areal);
friend areal operator/(areal, double);
/* mgn=luvun suuruus */
friend double mgn(areal c);
friend void mem00set(void*, int);
};
void mem00set(void *ptr, int kpl)
{
char *c=(char*)ptr;
for (int i=0; i<kpl; i++)
c[i]=(char)0;
}
int strToInt(char *m)
{
int x[8]; // max=2^16
int sum=0, i, j, k=0x01;
mem00set(x, sizeof(int)*8);
for (i=0, j=0; m[i]; i++)
if (m[i]>='0'&&m[i]<='9')
{
x[j]=(int)(m[i]-'0');
k*=(int)10;
++j;
}
for (i=0, k/=10; k; i++, k/=10)
sum+=(int)(x[i]*k);
return sum;
}
areal::areal(void)
{
mem00set(this, sizeof(areal));
}
areal::areal(double *a)
{
for (int i=0; i<DIM; i++)
f[i]=a[i];
}
////////////////////////////////////////////////////////////
// Muodostin areal(char*, ...) //
// //
// Esim. //
// //
// areal x("2d", 2, 3), x saa arvon (2.00000, 3.00000) //
// areal y("2f", 1.23, 4.56), y saa arvon (1.23, 4.56) //
// //
// - formaatissa annetaan ensin lukujen määrä, esim 4. //
// - d tarkoittaa, että luvut annetaan kokonaislukuina. //
// - f tarkoittaa, että luvut annetaan liukulukumuodossa. //
////////////////////////////////////////////////////////////
areal::areal(char *m, ...)
{
int sizeofint=sizeof(int);
int i, j, kpl=strToInt(m), step;
mem00set((void*)this, sizeof(areal));
unsigned long P=(unsigned long)&m+sizeof(char*);
for (j=0; m[j]; j++)
{
if (m[j]=='d') step=sizeof(int);
if (m[j]=='f') step=sizeof(double);
}
step=step==sizeofint? sizeofint: sizeof(double);
for (i=j=0x00; i<kpl; i++, P+=step)
f[i] = step == sizeof(int)?
*(int*)P: *(double*)P;
}
double abs(double x)
{
return x<0.0? -x: x;
}
double abs(areal a)
{
double sum=0.00;
for (int i=0; i<DIM; i++)
{
if (i & 1) sum-=a.f[i];
else /*~*/ sum+=a.f[i];
}
return abs(sum);
}
double mgn(areal a)
{
double sum=0.00;
for (int i=0; i<DIM; i++)
sum+=abs(a.f[i]);
return sum;
}
areal operator-(areal a)
{
for (int i=0; i<DIM; i++)
a.f[i]=-a.f[i];
return a;
}
areal operator+(areal a, areal b)
{
for (int i=0; i<DIM; i++)
a.f[i]+=b.f[i];
return a;
}
areal operator-(areal a, areal b)
{
for (int i=0; i<DIM; i++)
a.f[i]-=b.f[i];
return a;
}
areal operator*(areal a, areal b)
{
double x[DIM*2];
double y[DIM*2];
double t[DIM*2];
mem00set(x, sizeof(double)*DIM*2);
mem00set(y, sizeof(double)*DIM*2);
mem00set(t, sizeof(double)*DIM*2);
for (int i=0; i<DIM; i++)
{
x[i*2]=a.f[i];
y[i*2]=b.f[i];
}
for (int i=0; i<DIM*2; i++)
for (int j=0; j<DIM*2; j++)
t[i^j]+=x[i]*y[j];
for (int i=0; i<DIM; i++)
{
t[i]=t[i*2]-t[i*2+1];
}
return areal(t);
}
areal operator*(double k, areal b)
{
// muodostin antaa a:lle arvon k
areal a("1f", k);
return a*b;
}
areal operator*(areal a, double k)
{
// muodostin antaa b:lle arvon k
areal b("1f", k);
return a*b;
}
areal operator/(double k, areal b)
{
// muodostin antaa a:lle arvon k
areal a("1f", k);
return a/b;
}
areal operator/(areal a, double k)
{
// muodostin antaa b:lle arvon k
areal b("1f", k);
return a/b;
}
areal operator/(areal u, areal v)
{
double x[DIM*2];
double y[DIM*2];
mem00set(x, sizeof(double)*DIM*2);
mem00set(y, sizeof(double)*DIM*2);
for (int i=0; i<DIM; i++)
{
x[i*2]=u.f[i];
y[i*2]=v.f[i];
}
double X[DIM*2];
#define DIM2 (DIM*2)
double K[DIM2*DIM2+DIM2];
for (int n=0, p=0; n<DIM*2; n++)
{
for (int i=0; i<DIM*2; i++)
for (int j=0; j<DIM*2; j++)
{
if (n == (i^j))
{
K[p]=y[j]; ++p;
}
}
K[p]=x[n]; ++p;
}
// gsolve ratkaisee lineaarisen yhtälöryhmän
// Gaussin eliminointimenetelmällä...
int gsolve(double*, double*, int);
gsolve(X, (double*)K, DIM*2);
for (int i=0; i<DIM; i++)
X[i]=X[i*2]-X[i*2+1];
return areal(X);
}
////////////////////////////////////////////////
// qsolve palautusarvo: //
// 1 <=> yhtälöryhmälle oli ratkaisu //
// 0 <=> yhtälöryhmälle ei ollut ratkaisua //
////////////////////////////////////////////////
int gsolve(double *X, double *f, int n)
{
double k, J;
double *max, *mux;
int i, j, x, M=(int)n+1;
void swapRows(double*, double*, int);
for (i=0; i<n-1; i++)
{
max=f+i+i*M;
for (j=i+0x01; j < n; j++)
if (abs(*(mux=f+i+j*M))>abs(*max)) max=mux;
if (max!=(mux=f+i+i*M))
swapRows(max, mux, M-i);
if (abs(J=f[i+i*M])<(double)1e-15)
return 0; /* eipä ole ratkaisuva */
for (j=i+1; j<n; j++)
{
if (f[i+j*M])
{
k=-f[i+j*M]/J;
for (x=i+1; x<M; x++)
f[x+j*M]=f[x+j*M]+k*f[x+i*M];
}
}
}
for (i=n-1; i>=0; i--)
{
k=(double)f[n+i*M];
for (j=n-1; j > i; j--)
k=k-(double)(X[j]*f[j+i*M]);
if (abs(J=f[i+i*M])<(double)1e-15)
return 0; /* eipä ole ratkaisua */
else X[i]=k/J; /* onpa ratkaisu */
}
return 1;
}
void swapRows(double *a, double *b, int n)
{
for (int i=0; i<n; i++)
{
double tmp=a[i];
a[i] = b[i];
b[i]=tmp;
}
}
void print(areal a)
{
for (int i=0; i<DIM; i++)
{
if (!i) printf("(");
printf("%0.5f", a.f[i]);
if (i+1==DIM) printf(")\n");
else printf(", ");
}
}
/* ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** */
/* ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** */
/* ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** DEBUG ** *** *** ** */
// tuotetaan satunnaisluku
unsigned long rnd(unsigned long max)
{
static unsigned long R16A=1L;
static unsigned long R16B=2L;
R16A -= 0x012357bfL;
R16B += 0xfedca201L;
R16A += (R16A>>16)^(R16B<<16);
R16B -= (R16A<<16)^(R16B>>16);
return (R16A^R16B)%max;
}
// tuotetaan satunnaisluku
double drand(void)
{
unsigned long jak=0xfffffffb;
double sg=rnd(0x02)? -1: 1;
double x=(double)rnd(jak);
return sg * x/(double)jak;
}
areal areal_rnd(void)
{
double x[DIM];
for (int i=0; i<DIM; i++)
x[i]=drand();
return areal(x);
}
#pragma argsused
int main(int kpl, char* asia[])
{
for (int i=0; i<16; i++)
{
areal a=areal_rnd();
areal b=areal_rnd();
printf("ehto 1: a*b = b*a \n");
print(a*b);
print(b*a);
printf("\n");
printf("ehto 2: a/b=c <=> a=b*c \n");
areal c=a/b;
print(a);
print(b*c);
printf("\n");
printf("ehto 3: abs(a)*abs(b) = abs(a*b) \n");
printf("%0.9f\n", abs(a)*abs(b));
printf("%0.9f\n", abs(a*b));
printf("\n");
}
return 0;
}
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
|
Läskiperse
|
Lähetetty: Pe Heinä 13, 2012 3:50 pm |
|
Liittynyt: La Joulu 11, 2010 6:24 pm
Viestit: 491
|
|
tuota...ei sillä varmaan ole väliä, mutta koodaan copy-pastella. jokaisessa uudessa painoksessa pyrin tiivistämään ja selkeyttämään koodia. debugaan joka funktion erikseen, ja lopuksi yhdessä.
olen sitä mieltä, että koodin pitää kertoa tiiviinä algoritmina, mitä tapahtuu. varsinkin matemaattisissa koodeissa ei juuri kommentointia kaivata, koska se vain häiritsee koodin lukemista. funktio ei saa olla juuri A4:sta pidempi, koska silloin tulee mössöä - eikä edes runsas kommentointi paranna asiaa.
NMP:llä koodit olivat niin p:stä, että otin itse lopputilin. muissa työpaikoissa on ollut koodien suhteen sama meininki. työpaikkojen merkitys on siinä, että pidetään runsaat kesäjuhlat ja pikkujoulut. muuten haistatan tietotekniikalle (ATK) pitkät p:t.
pitää seuraavaksi liittää polynomi areal:iin, sekä tehdä mainiin syy-seuraus:
syy0=(a0, b0, c0, 0, 0, 0, 0, 0), seuraus0=(d0, e0, f0, g0, h0, i0, j0, k0)
syy1=(a1, b1, c1, 0, 0, 0, 0, 0), seuraus0=(d1, e1, f1, g1, h1, i1, j1, k1)
syy2=(a2, b2, c2, 0, 0, 0, 0, 0), seuraus0=(d2, e2, f2, g2, h2, i2, j2, k2)
...
lisään kohinan randomilla, ja lasken käyräparven lyhyeksi polynomiksi. eräs keskustelija virkkoi, etten oikein ymmärrä kohinaa. minun mielestä cern:ssä ei ymmäretä kohinaa sen enempää. yksi työkaveri oli ollut cern:ssä töissä. ihmettelin, kun kaveri ei osannut laskea. tai ei sen väliä, pääasia, että palkka juoksee.
_________________
...haaveet suomalaisesta Läskiperseestä matematiikan historian aikakirjoissa saa kyllä haudata. -Cargo
|
|
| Ylös |
|
|
Paikallaolijat |
Käyttäjiä lukemassa tätä aluetta: visti ja 11 vierailijaa |
|
Et voi kirjoittaa uusia viestejä
Et voi vastata viestiketjuihin
Et voi muokata omia viestejäsi
Et voi poistaa omia viestejäsi
|
|
|
2
Kirjaudu sisään
Rekisteröidy
UKK
Etsi
