procedure Sosu(var intArr :array of Integer); var sosuflag:array[word]of boolean; //素数かどうかを判定するフラグ i,j :Integer; begin fillchar(sosuflag,SizeOf(Sosuflag),1); //まずすべての数を素数と見なしてから後で、割り切れるものをはじく。 i := 2; //素数の最小値 repeat j := i ; inc(j,i); while j <= $FFFF do begin sosuflag[j] := false; //iの倍数なので素数ではない Inc(J,I); end; repeat Inc(i); until sosuflag[i]; //次に小さい素数を探す until i > $8FFF; j := 1; for i := 0 to High(IntArr) do begin repeat inc(j); if j > $FFFF then raise Exception.Create('大きすぎて素数が求められません'); until sosuflag[j]; IntArr[i] := J; end; end; 0819１３２人目の素数さん2020/11/28(土) 08:59:49.93ID:LpYp+oBb 悪霊退散！！！
LET n = 100 DIM s(n) MAT s = ZER ! 配列 s の全要素に 0(zero) を代入 ※s = 0 ではダメ LET k = 0 FOR i = 2 TO n IF s(i) = 0 THEN PRINT USING "####":i; LET k = k + 1 IF MOD(k,10) = 0 THEN PRINT END IF FOR j = i^2 TO N STEP i LET s(j) = 1 NEXT j END IF NEXT i 0822日高2020/11/28(土) 09:04:47.02ID:0fpuH75L (修正10) 【定理】n>2のとき、x^n+y^n=z^nのx,y,zは自然数とならない。 【証明】x^n+y^n=z^nを、z=x+rとおいてx^n+y^n=(x+r)^n…(1)とする。 (1)をr^(n-1){(y/r)^n-1}=an{x^(n-1)+…+r^(n-2)x}(1/a)…(2)と変形する。 (2)はa=1、r^(n-1)=nのとき、x^n+y^n=(x+n^{1/(n-1)})^n…(3)となる。 (2)はa=1以外、r^(n-1)=anのとき、x^n+y^n=(x+(an)^{1/(n-1)})^n…(4)となる。 (3)はn=2のとき、x,y,zは整数比となりえる。n>2のとき、x,y,zは整数比とならない。 (4)のx,y,zは(3)のx,y,zのa^{1/(n-1)}倍となる。 ∴n>2のとき、x^n+y^n=z^nのx,y,zは自然数とならない。 0823１３２人目の素数さん2020/11/28(土) 09:10:01.88ID:LpYp+oBb 悪霊退散！！！
REM 既約なピタゴラス数 FUNCTION gcd(a,b) DO WHILE b <> 0 LET r = MOD(a,b) LET a = b LET b = r LOOP LET gcd = a END FUNCTION
LET LAST = 200
REM ピタゴラス数を求める FOR x = 1 TO LAST FOR y = x + 1 TO LAST LET z = SQR(x^2+y^2) IF INT(z) = z THEN IF gcd(x,y) = 1 AND gcd(x,z) = 1 AND gcd(y,z) = 1 THEN PRINT USING "##### ##### #####": x,y,z END IF END IF NEXT y NEXT x
END 0824日高2020/11/28(土) 09:11:49.19ID:0fpuH75L >800 x,y,zに関する方程式(3)の満たすべき条件は 「x^n +y^n=z^n」だけではない。 という意味です。
どの部分が、フニャフニャ回避でしょうか？ 0828１３２人目の素数さん2020/11/28(土) 09:21:50.04ID:LpYp+oBb 悪霊退散！！！ begin with PaintScreen1 do begin OffBmp.Canvas.Pen.Width := SubWidth; OffBmp.Canvas.Pen.Color := TColor($971D4F);
//初期値 h := StrToFloat(Edh.Text); Et := StrToFloat(Edt0.Text); Ev := StrToFloat(Edv0.Text); Ex := StrToFloat(Edx0.Text); OffBmp.Canvas.MoveTo( RealToDispX(Et),RealToDispY(Ex) ); for i := 0 to 100 do begin 0829１３２人目の素数さん2020/11/28(土) 09:22:24.03ID:LpYp+oBb 悪霊退散！！！ //ルンゲ・クッタ k1 := h*DFunc1(Et, Ex, Ev); m1 := h*DFunc2(Et, Ex, Ev);
どの部分が、自分の論理が破綻してしまうような質問でしょうか？ 0831１３２人目の素数さん2020/11/28(土) 09:23:01.09ID:LpYp+oBb 悪霊退散！！！ //初期値 h := StrToFloat(Edh.Text); Et := StrToFloat(Edt0.Text); Ev := StrToFloat(Edv0.Text); Ex := StrToFloat(Edx0.Text); OffBmp.Canvas.MoveTo( RealToDispX(Et),RealToDispY(Ex) ); for i := 0 to 100 do //負方向の計算 begin k1 := -h*DFunc1(Et, Ex, Ev); m1 := -h*DFunc2(Et, Ex, Ev);
どの部分が、己の論理の穴をちゃんと認識していることに、なるのでしょうか？ 0833１３２人目の素数さん2020/11/28(土) 09:26:25.49ID:LpYp+oBb 悪霊退散！！！ procedure fft; var K,L,KD : Integer ; wc,ws : Array of Extended ; procedure fftint ; var s : Integer ; wk : Extended ; begin For s := 0 to KD-1 do begin wk := 2.0 * Pi * s / K ; wc[s] := Cos(wk) ; ws[s] := -Sin(wk) ; end ; end ; function bitrev(ip : Integer) : Integer ; var i,w : Integer ; begin w := 0 ; For I := 1 to L do begin w := w * 2 + (ip mod 2) ; ip := ip div 2 ; end ; Result := w ; end ; procedure cfft(inv : Integer) ; var i,j,Li,sn,i0,i1,expon,iw : Integer ; wk,yr,yi,sign,wwc,wws : Extended ; 0834１３２人目の素数さん2020/11/28(土) 09:27:01.68ID:LpYp+oBb 悪霊退散！！！ begin if inv = 1 then sign := 1.0 else sign := -1.0 ; Li := K ; iw := 1 ; For I := 1 to L do begin iw := iw * 2 ; Li := Li div 2 ; sn := 0 ; while sn < K-1 do begin For J := 0 to Li-1 do begin expon := (bitrev(sn) mod iw) * Li ; wws := sign * ws[expon] ; wwc := wc[expon] ; i0 := sn ; i1 := i0 + Li ; yr := xr[i1] * wwc - xi[i1] * wws ; yi := xr[i1] * wws + xi[i1] * wwc ; xr[i1] := xr[i0] - yr ; xi[i1] := xi[i0] - yi ; xr[i0] := xr[i0] + yr ; xi[i0] := xi[i0] + yi ; Inc(sn) ; end ; Inc(sn,Li) ; end ; end ; 0835１３２人目の素数さん2020/11/28(土) 09:27:27.58ID:LpYp+oBb 悪霊退散！！！ For i := 0 to K-1 do begin j := bitrev(i) ; if inv = 1 then begin cr[j] := xr[i] / K ; ci[j] := xi[i] / K ; end else begin cr[j] := xr[i] ; ci[j] := xi[i] ; end ; end ; end ; var s : Integer ; begin //分割数を得る K := 1024 ; //2^nの形 L := Trunc(log2(K)) ; KD := K div 2 ; //データはxr,xiで渡される。 SetLength(wc,KD) ; SetLength(ws,KD) ; //手続き呼び出し fftint ; //以下で逆高速フーリエ変換ならば引数を-1にすれば大丈夫 cfft(1) ; end; 0836１３２人目の素数さん2020/11/28(土) 09:29:56.61ID:LpYp+oBb 悪霊退散！！！ //ガウス・ジョルダン〜Pivot選択あり：逆行列ルーティン付き function GaussJordanPv(N: Integer):Integer; var pRow,pv, k, j: Integer; mMax,R_pivot, temp: Extended; begin //単位行列の設定 for k := 1 to N do for j := 1 to N do if k = j then RA[k][j] := 1.0 else RA[k][j] := 0.0; for pv := 1 to N do //行ループ begin mMax := 0.000001; for k := pv to N do //行ループ 最大値探索 begin if Abs(A[k][pv]) > mMax then begin mMax := Abs(A[k][pv]); pRow := k; end; end; if mMax <= 0.000001 then //誤差対策 begin MessageDlg('解が存在しないかまたは不定です！', mtwarning, [mbok], 0); Result := 0; Exit; end; 0837１３２人目の素数さん2020/11/28(土) 09:30:39.38ID:LpYp+oBb 悪霊退散！！！ //行の入れ替え if pv <> pRow then begin for k := 1 to N+1 do //列ループ begin temp := A[pv][k]; A[pv][k] := A[pRow][k]; A[pRow][k] := temp; end; for k := 1 to N do //列ループ 単位行列 begin temp := RA[pv][k]; RA[pv][k] := RA[pRow][k]; RA[pRow][k] := temp; end; end;
//ピボット行の処理 ⇒ 対角成分 = 1 R_pivot := 1.0/A[pv][pv];//ピボットの逆数 for j := 1 to N+1 do //列ループ A[pv][j] := A[pv][j]*R_pivot; for j := 1 to N do //列ループ 単位行列 RA[pv][j] := RA[pv][j]*R_pivot; //ピボット行以外の処理 ⇒ ピボット列 = 0 for k := 1 to N do begin temp := A[k][pv]; //消去する係数 begin for j := pv to N+1 do //ピボット列以降を処理 if k <> pv then A[k][j] := A[k][j] - temp*A[pv][j]; for j := 1 to N do //全列処理（単位行列） if k <> pv then RA[k][j] := RA[k][j] - temp*RA[pv][j]; end; end; end; Result := 1; end; 0838日高2020/11/28(土) 09:30:40.59ID:0fpuH75L >811 nが自然数でなければ明確に反例が知られている。
nが無理数ならば、反例があります。 0839１３２人目の素数さん2020/11/28(土) 09:33:10.91ID:LpYp+oBb 悪霊退散！！！ function Formalize(var S: String): Boolean; var i : Integer; dmy: String; begin Result := False; if S = '' then Exit; // 空文字はエラー for i := 1 to Length(S) do // 不正な文字があるとエラー if not ( S[i] in ['0'..'9', '+', '-', '*', '/', '(', ')', '.'] ) then Exit; dmy := S; i := 1; while i <= Length(S) do if not ( S[i] in ['(', ')'] ) then Delete(S, i, 1) else Inc(i); while Pos('()', S) > 0 do Delete(S, Pos('()', S), 2); if Length(S) > 0 then Exit; S := dmy; if S[1] in ['+', '-', '*', '/', ')', '.'] then Exit; if S[Length(S)] in ['+', '-', '*', '/', '(', '.'] then Exit; for i := 1 to Length(S) - 1 do begin if (S[i] in ['+', '-', '*', '/', '.', '(']) and (S[i+1] in ['+', '-', '*', '/', '.', ')']) then Exit; end; for i := 2 to Length(S) - 1 do if (S[i] = '.') then if not ((S[i-1] in ['0'..'9']) and (S[i+1] in ['0'..'9'])) then Exit; Result := True; end; 0840１３２人目の素数さん2020/11/28(土) 09:34:08.30ID:LpYp+oBb 悪霊退散！！！ procedure NextToken; begin case GNum of '0'..'9': GetNumber; GTOKEN := C_NUMBER; '+','-': GOP := GNum; GTOKEN := C_ADD; '*','/': GOP := GNum; GTOKEN := C_MUL; '(': GTOKEN := C_LPAREN; ')': GTOKEN := C_RPAREN; ')': GTOKEN := C_RPAREN; '#': GTOKEN := C_OTHERS; end; if not (GNum in ['0'..'9', '#']) then GNum := ReadChar; // 先読み end; 0841日高2020/11/28(土) 09:34:13.71ID:0fpuH75L (修正10) 【定理】n>2のとき、x^n+y^n=z^nのx,y,zは自然数とならない。 【証明】x^n+y^n=z^nを、z=x+rとおいてx^n+y^n=(x+r)^n…(1)とする。 (1)をr^(n-1){(y/r)^n-1}=an{x^(n-1)+…+r^(n-2)x}(1/a)…(2)と変形する。 (2)はa=1、r^(n-1)=nのとき、x^n+y^n=(x+n^{1/(n-1)})^n…(3)となる。 (2)はa=1以外、r^(n-1)=anのとき、x^n+y^n=(x+(an)^{1/(n-1)})^n…(4)となる。 (3)はn=2のとき、x,y,zは整数比となりえる。n>2のとき、x,y,zは整数比とならない。 (4)のx,y,zは(3)のx,y,zのa^{1/(n-1)}倍となる。 ∴n>2のとき、x^n+y^n=z^nのx,y,zは自然数とならない。 0842１３２人目の素数さん2020/11/28(土) 09:35:22.57ID:LpYp+oBb 悪霊退散！！！ function Expression: Extended; var op : Char sign: Integer; u,v : Extended; begin u := Term; while GTOKEN := C_ADD do begin op := GOP; // オペレータを保存 NextToken; // '+','-' を読みとばす v := Term; if op = '+' then u := u+v else u := u-v; end; Result := u; end;
function Term: Extended; var op: Char; u,v: Extended; begin u := Factor;; while GTOKEN := C_MUL do begin op := GOP; // オペレータを保存 NextToken; // '*','/' を読みとばす v := Factor; case op of '*': u := u*v; '/': u := u/ v; else //Error! end; end; Result := u; end; 0843１３２人目の素数さん2020/11/28(土) 09:36:29.99ID:LpYp+oBb 悪霊退散！！！ function Factor: Extended; var v: Extended; begin case GTOKEN of C_LPAREN: // 左括弧の場合 NextToken; // '(' を読みとばす v := Expression; // 「式」の処理 if GTOKEN = C_RPAREN then // ')' が来ているはず。チェック NextToken // ')' を読みとばす else ErrorOut; C_NUMBER: // 数値の場合 v := GVALUE; // 数値の処理をする NextToken; // 数値を読みとばす else: ErrorOut; // "(" でも数値でもなければ、エラー end; Result := v; end; 0844１３２人目の素数さん2020/11/28(土) 09:39:14.00ID:LpYp+oBb 悪霊退散！！！ int main() { int i; double k, p, dk, kmin, kmax, pmin, pmax;
printf("kmin = "); scanf("%lf", &kmin); printf("kmax = "); scanf("%lf", &kmax); printf("pmin = "); scanf("%lf", &pmin); printf("pmax = "); scanf("%lf", &pmax); gr_on(); gr_window(kmin, pmin, kmax, pmax, 0, 0); dk = (kmax - kmin) / (XMAX - 1); for (k = kmin; k <= kmax; k += dk) { p = 0.3; for (i = 1; i <= 50; i++) p += k * p * (1 - p); for (i = 51; i <= 100; i++) { if (p >= pmin && p <= pmax) gr_wdot(k, p, WHITE); p += k * p * (1 - p); } } hitanykey(); return EXIT_SUCCESS; } 0845１３２人目の素数さん2020/11/28(土) 09:40:20.75ID:LpYp+oBb 悪霊退散！！！ #define N 24 #define PI 3.14159265358979323846264 long double latan(long double x) /* アークタンジェント */ { int i, sgn; long double a;
if (x > 1) { sgn = 1; x = 1 / x; } else if (x < -1) { sgn = -1; x = 1 / x; } else sgn = 0; a = 0; for (i = N; i >= 1; i--) a = (i * i * x * x) / (2 * i + 1 + a); if (sgn > 0) return PI / 2 - x / (1 + a); if (sgn < 0) return -PI / 2 - x / (1 + a); /* else */ return x / (1 + a); }
void corrcoef1(int n, float x[], float y[]) { int i; float sx, sy, sxx, syy, sxy, dx, dy;
sx = sy = sxx = syy = sxy = 0; for (i = 0; i < n; i++) { sx += x[i]; sy += y[i]; } sx /= n; sy /= n; for (i = 0; i < n; i++) { dx = x[i] - sx; dy = y[i] - sy; sxx += dx * dx; syy += dy * dy; sxy += dx * dy; } sxx = sqrt(sxx / (n - 1)); syy = sqrt(syy / (n - 1)); sxy /= (n - 1) * sxx * syy; printf("標準偏差 %g %g 相関係数 %g\n", sxx, syy, sxy); }
void corrcoef2(int n, float x[], float y[]) { int i; float sx, sy, sxx, syy, sxy;
sx = sy = sxx = syy = sxy = 0; for (i = 0; i < n; i++) { sx += x[i]; sy += y[i]; sxx += x[i] * x[i]; syy += y[i] * y[i]; sxy += x[i] * y[i]; } sx /= n; sxx = (sxx - n * sx * sx) / (n - 1); sy /= n; syy = (syy - n * sy * sy) / (n - 1); if (sxx > 0) sxx = sqrt(sxx); else sxx = 0; if (syy > 0) syy = sqrt(syy); else syy = 0; sxy = (sxy - n * sx * sy) / ((n - 1) * sxx * syy); printf("標準偏差 %g %g 相関係数 %g\n", sxx, syy, sxy); } 0850１３２人目の素数さん2020/11/28(土) 09:46:21.82ID:LpYp+oBb 悪霊退散！！！
void corrcoef3(int n, float x[], float y[]) { int i; float sx, sy, sxx, syy, sxy, dx, dy;
sx = sy = sxx = syy = sxy = 0; for (i = 0; i < n; i++) { dx = x[i] - sx; sx += dx / (i + 1); dy = y[i] - sy; sy += dy / (i + 1); sxx += i * dx * dx / (i + 1); syy += i * dy * dy / (i + 1); sxy += i * dx * dy / (i + 1); } sxx = sqrt(sxx / (n - 1)); syy = sqrt(syy / (n - 1)); sxy /= (n - 1) * sxx * syy; printf("標準偏差 %g %g 相関係数 %g\n", sxx, syy, sxy); } 0851１３２人目の素数さん2020/11/28(土) 09:48:42.35ID:LpYp+oBb 悪霊退散！！！ #include <math.h> #define N 8 static double coef[20] = { 8.333333333333333333333333333e-2, /* 1/12 */ -1.388888888888888888888888889e-3, /* -1/720 */ 3.306878306878306878306878307e-5, /* 1/30240 */ -8.267195767195767195767195767e-7, /* -1/1209600 */ 2.087675698786809897921009032e-8, /* 1/47900160 */ -5.284190138687493184847682202e-10, 1.338253653068467883282698098e-11, -3.389680296322582866830195391e-13, 8.586062056277844564135905450e-15, -2.174868698558061873041516424e-16, 5.509002828360229515202652609e-18, -1.395446468581252334070768626e-19, 3.534707039629467471693229977e-21, -8.953517427037546850402611251e-23, 2.267952452337683060310950058e-24, -5.744790668872202445263829503e-26, 1.455172475614864901866244572e-27, -3.685994940665310178130050728e-29, 9.336734257095044668660153106e-31, -2.365022415700629886484029550e-32 }; 0852１３２人目の素数さん2020/11/28(土) 09:49:07.89ID:LpYp+oBb 悪霊退散！！！ double zeta(double x) { int i; double powNx, w, z, zprev;
z = 1; for (i = 2; i < N; i++) { zprev = z; z += pow(i, -x); if (z == zprev) return z; } powNx = pow(N, x); w = x / (N * powNx); z += 0.5 / powNx + N / ((x - 1) * powNx) + coef[0] * w; for (i = 1; i < 20 && z != zprev; i++) { w *= (x + 2 * i - 1) * (x + 2 * i) / (N * N); zprev = z; z += coef[i] * w; } return z; } 0853１３２人目の素数さん2020/11/28(土) 09:50:28.37ID:LpYp+oBb 悪霊退散！！！ unsigned phi(unsigned x) { unsigned d, t;
t = x; if (x % 2 == 0) { t /= 2; do { x /= 2; } while (x % 2 == 0); } d = 3; while (x / d >= d) { if (x % d == 0) { t = t / d * (d - 1); do { x /= d; } while (x % d == 0); } d += 2; } if (x > 1) t = t / x * (x - 1); return t; }
i = 0; j = N - 1; while (i < j) { k = (i + j) / 2; if (x[k] < t) i = k + 1; else j = k; } if (i > 0) i--; h = x[i + 1] - x[i]; d = t - x[i]; return (((z[i + 1] - z[i]) * d / h + z[i] * 3) * d + ((y[i + 1] - y[i]) / h - (z[i] * 2 + z[i + 1]) * h)) * d + y[i]; } 0856１３２人目の素数さん2020/11/28(土) 09:53:24.39ID:LpYp+oBb 悪霊退散！！！ #define TEST 1
#if TEST int count = 0; #endif
int A(int x, int y) { #if TEST count++; #endif if (x == 0) return y + 1; if (y == 0) return A(x - 1, 1); return A(x - 1, A(x, y - 1)); }
void encode(void) /* 圧縮 */ { int c; unsigned long range, maxcount, incount, cr, d; unsigned short low, high; static unsigned long count[N];
for (c = 0; c < N; c++) count[c] = 0; /* 頻度の初期化 */ while ((c = getc(infile)) != EOF) count[c]++; /* 各文字の頻度 */ incount = 0; maxcount = 0; /* 原文の大きさ, 頻度の最大値 */ for (c = 0; c < N; c++) { incount += count[c]; if (count[c] > maxcount) maxcount = count[c]; } if (incount == 0) return; /* 0バイトのファイル */ /* 頻度合計が {\tt Q1} 未満, 各頻度が1バイトに収まるよう規格化 */ d = max((maxcount + N - 2) / (N - 1), (incount + Q1 - 257) / (Q1 - 256)); if (d != 1) for (c = 0; c < N; c++) count[c] = (count[c] + d - 1) / d; cum[0] = 0; for (c = 0; c < N; c++) { fputc((int)count[c], outfile); /* 頻度表の出力 */ cum[c + 1] = cum[c] + (unsigned)count[c]; /* 累積頻度 */ } 0858１３２人目の素数さん2020/11/28(土) 09:56:31.59ID:LpYp+oBb 悪霊退散！！！ outcount = N; rewind(infile); incount = 0; /* 巻き戻して再走査 */ low = 0; high = USHRT_MAX; ns = 0; while ((c = getc(infile)) != EOF) { /* 各文字を符号化 */ range = (unsigned long)(high - low) + 1; high = (unsigned short) (low + (range * cum[c + 1]) / cum[N] - 1); low = (unsigned short) (low + (range * cum[c ]) / cum[N]); for ( ; ; ) { if (high < Q2) output(0); else if (low >= Q2) output(1); else if (low >= Q1 && high < Q3) { ns++; low -= Q1; high -= Q1; } else break; low <<= 1; high = (high << 1) + 1; } if ((++incount & 1023) == 0) printf("%12lu\r", incount); } 0859１３２人目の素数さん2020/11/28(土) 09:57:10.88ID:LpYp+oBb 悪霊退散！！！ int binarysearch(unsigned x) /* $\mbox{\tt cum[i]} \le x < \mbox{\tt cum[i+1]}$ となる {\tt i} を二分探索で求める */ { int i, j, k;
i = 1; j = N; while (i < j) { k = (i + j) / 2; if (cum[k] <= x) i = k + 1; else j = k; } return i - 1; } 0860１３２人目の素数さん2020/11/28(土) 09:57:37.15ID:LpYp+oBb 悪霊退散！！！ void decode(unsigned long size) /* 復元 */ { int c; unsigned char count[N]; unsigned short low, high, value; unsigned long i, range;
if (size == 0) return; /* 0バイトのファイル */ cum[0] = 0; for (c = 0; c < N; c++) { count[c] = fgetc(infile); /* 頻度分布を読む */ cum[c + 1] = cum[c] + count[c]; /* 累積頻度を求める */ } value = 0; for (c = 0; c < USHRT_BIT; c++) value = 2 * value + getbit(); /* バッファを満たす */ low = 0; high = USHRT_MAX; for (i = 0; i < size; i++) { /* 各文字を復元する */ range = (unsigned long)(high - low) + 1; c = binarysearch((unsigned)((((unsigned long) (value - low) + 1) * cum[N] - 1) / range)); high = (unsigned short) (low + (range * cum[c + 1]) / cum[N] - 1); low = (unsigned short) (low + (range * cum[c ]) / cum[N]); for ( ; ; ) { if (high < Q2) { /* 何もしない */ } else if (low >= Q2) { /* 何もしない */ } else if (low >= Q1 && high < Q3) { value -= Q1; low -= Q1; high -= Q1; } else break; low <<= 1; high = (high << 1) + 1; value = (value << 1) + getbit(); /* 1ビット読む */ } putc(c, outfile); /* 復元した文字を書き出す */ if ((i & 1023) == 0) printf("%12lu\r", i); } printf("%12lu\n", size); /* 原文のバイト数 */ } 0861１３２人目の素数さん2020/11/28(土) 10:07:04.52ID:LpYp+oBb 悪霊退散！！！ int change(int n, int k) /* 再帰版 */ { int s;
if (n < 0) return 0; s = 1 + n / 5 + change(n - 10, 10); if (k >= 50) s += change(n - 50, 50); if (k >= 100) s += change(n - 100, 100); return s; }
int change1(int n) /* 非再帰版 */ { int i, j, s, t, u;
s = 0; for (i = n / 100; i >= 0; i--) { /* 100円玉 */ t = n - 100 * i; for (j = t / 50; j >= 0; j--) { /* 50円玉 */ u = t - 50 * j; s += (1 + u / 5 - u / 10) * (1 + u / 10); } } return s; }
#include <stdio.h> #include <stdlib.h>
int main() { int i;
printf("お金の払い方\n"); printf(" 金額 再帰版 非再帰版\n"); for (i = 0; i <= 500; i += 5) printf("%6d %8d %8d\n", i, change(i, i), change1(i)); return EXIT_SUCCESS; } 0862１３２人目の素数さん2020/11/28(土) 10:08:41.83ID:LpYp+oBb #include <math.h> #define PI 3.14159265358979323846264 double p_nor(double z) /* 正規分布の下側累積確率 */ { int i; double z2, prev, p, t;
z2 = z * z; t = p = z * exp(-0.5 * z2) / sqrt(2 * PI); for (i = 3; i < 200; i += 2) { prev = p; t *= z2 / i; p += t; if (p == prev) return 0.5 + p; } return (z > 0); }
s = EULER + log(x); x = - x * x; t = 1; for (k = 2; k < 1000; k += 2) { t *= x / ((k - 1) * k); u = s; s += t / k; if (s == u) return s; } printf("Si_series(): 収束しません.\n"); return s; }
printf("n = "); scanf("%lu", &n); while (n > 1) { if (n & 1) { /* 奇数 */ if (n > LIMIT) { printf("\nOverflow\n"); return 1; } else n = 3 * n + 1; } else n /= 2; printf(" %lu", n); } printf("\n"); return EXIT_SUCCESS; } 0865１３２人目の素数さん2020/11/28(土) 10:11:57.73ID:LpYp+oBb 悪霊退散！！！ int comb(int n, int k) { if (k == 0 || k == n) return 1; /* if (k == 1) return n; */ return comb(n - 1, k - 1) + comb(n - 1, k); }
unsigned long combination(int n, int k) { int i, j; unsigned long a[17];
if (n - k < k) k = n - k; if (k == 0) return 1; if (k == 1) return n; if (k > 17) return 0; /* error */ for (i = 1; i < k; i++) a[i] = i + 2; for (i = 3; i <= n - k + 1; i++) { a[0] = i; for (j = 1; j < k; j++) a[j] += a[j - 1]; } return a[k - 1]; } 0866１３２人目の素数さん2020/11/28(土) 10:13:43.71ID:LpYp+oBb #define SCALAR double #include "matutil.c" #include <math.h>
double lu(int n, matrix a, int *ip) /* LU分解 */ { int i, j, k, ii, ik; double t, u, det; vector weight;
weight = new_vector(n); /* {\tt weight[0..n-1]} の記憶領域確保 */ det = 0; /* 行列式 */ for (k = 0; k < n; k++) { /* 各行について */ ip[k] = k; /* 行交換情報の初期値 */ u = 0; /* その行の絶対値最大の要素を求める */ for (j = 0; j < n; j++) { t = fabs(a[k][j]); if (t > u) u = t; } if (u == 0) goto EXIT; /* 0 なら行列はLU分解できない */ weight[k] = 1 / u; /* 最大絶対値の逆数 */ } det = 1; /* 行列式の初期値 */ for (k = 0; k < n; k++) { /* 各行について */ u = -1; for (i = k; i < n; i++) { /* より下の各行について */ ii = ip[i]; /* 重み×絶対値 が最大の行を見つける */ t = fabs(a[ii][k]) * weight[ii]; if (t > u) { u = t; j = i; } } ik = ip[j]; if (j != k) { ip[j] = ip[k]; ip[k] = ik; /* 行番号を交換 */ det = -det; /* 行を交換すれば行列式の符号が変わる */ } u = a[ik][k]; det *= u; /* 対角成分 */ if (u == 0) goto EXIT; /* 0 なら行列はLU分解できない */ for (i = k + 1; i < n; i++) { /* Gauss消去法 */ ii = ip[i]; t = (a[ii][k] /= u); for (j = k + 1; j < n; j++) a[ii][j] -= t * a[ik][j]; } } EXIT: free_vector(weight); /* 記憶領域を解放 */ return det; /* 戻り値は行列式 */ } 0867１３２人目の素数さん2020/11/28(土) 10:14:24.10ID:LpYp+oBb 悪霊退散！！！ double matinv(int n, matrix a, matrix a_inv) { int i, j, k, ii; double t, det; int *ip; /* 行交換の情報 */
ip = malloc(sizeof(int) * n); if (ip == NULL) error("記憶領域不足"); det = lu(n, a, ip); if (det != 0) for (k = 0; k < n; k++) { for (i = 0; i < n; i++) { ii = ip[i]; t = (ii == k); for (j = 0; j < i; j++) t -= a[ii][j] * a_inv[j][k]; a_inv[i][k] = t; } for (i = n - 1; i >= 0; i--) { t = a_inv[i][k]; ii = ip[i]; for (j = i + 1; j < n; j++) t -= a[ii][j] * a_inv[j][k]; a_inv[i][k] = t / a[ii][i]; } } free(ip); return det; }
double infinity_norm(int n, matrix a) /* ∞ノルム */ { int i, j; double rowsum, max;
max = 0; for (i = 0; i < n; i++) { rowsum = 0; for (j = 0; j < n; j++) rowsum += fabs(a[i][j]); if (rowsum > max) max = rowsum; } return max; } 0868１３２人目の素数さん2020/11/28(土) 10:43:01.09ID:LpYp+oBb 悪霊退散!!! #include <stdio.h> #include <stdlib.h> #include <string.h> #include <ctype.h> enum {FALSE, TRUE};
#define N 10 /* 最大の行数 */
int imax, jmax, solution, word[N][128], digit[256], low[256], ok[10];
void found(void) /* 解の表示 */ { int i, j, c;
printf("\n解 %d\n", ++solution); for (i = 0; i <= imax; i++) { for (j = jmax; j >= 0; j--) { c = word[i][j]; if (c != '\0') printf("%d", digit[c]); else printf(" "); } printf("\n"); } } 0869１３２人目の素数さん2020/11/28(土) 10:43:37.17ID:LpYp+oBb 悪霊退散！！！
void try(int sum) /* 再帰的に試みる */ { static int i = 0, j = 0, carry; int c, d;
c = word[i][j]; if (i < imax) { i++; if ((d = digit[c]) < 0) { /* 定まっていないなら */ for (d = low[c]; d <= 9; d++) if (ok[d]) { digit[c] = d; ok[d] = FALSE; try(sum + d); ok[d] = TRUE; } digit[c] = -1; } else try(sum + d); i--; } else { j++; i = 0; d = sum % 10; carry = sum / 10; if (digit[c] == d) { if (j <= jmax) try(carry); else if (carry == 0) found(); } else if (digit[c] < 0 && ok[d] && d >= low[c]) { digit[c] = d; ok[d] = FALSE; if (j <= jmax) try(carry); else if (carry == 0) found(); digit[c] = -1; ok[d] = TRUE; } j--; i = imax; } } 0870１３２人目の素数さん2020/11/28(土) 10:44:57.61ID:LpYp+oBb 悪霊退散!!!
if (a > b) { t = a; a = b; b = t; } t = r * (b - a); c = a + t; d = b - t; fc = f(c); fd = f(d); for ( ; ; ) { if (fc > fd) { a = c; c = d; fc = fd; d = b - r * (b - a); if (d - c <= tolerance) return c; fd = f(d); } else { b = d; d = c; fd = fc; c = a + r * (b - a); if (d - c <= tolerance) return d; fc = f(c); } } }
#include <stdio.h> #include <stdlib.h> #define TEST 1
