くだすれFORTRAN（超初心者用）

くだすれFORTRAN（超� ..

39:デフォルトの名無しさん
06/01/28 20:02:47
とりあえずフォートランに多相型っていつはいるんですか？

40:デフォルトの名無しさん
06/01/28 20:10:27
>>39
polymorphismならFORTRAN2003でそれなりに入っている。

41:デフォルトの名無しさん
06/01/28 20:10:42
>>35
まず、両方の曲線の交点が２つあって、その間の面積は解析的に求めることは簡単ですね。
で、その解析的に出た面積はaが分かれば数値的に出てくると。
そのaは二元連立方程式の解の一つであると。

二元連立方程式はプログラミングしないでも、中学レベルの頭ならすぐ解ける。

どこにプログラミングの必要性があるのか疑問。

42:デフォルトの名無しさん
06/01/28 20:17:28
>>38
いまいち問題設定が見えてこないなｗｗｗ　これが問題のすべてか？背景は無いのか？

ガウスの消去法は、例外処理を考えなければ、手でやるとおりに組めｗ
例外処理はめんどくさい。

２＊２限定なら、高校で習った公式に突っ込んだほうが早い。例外がdet=0で排除できるし。
というか手でといたほうが早いわけだが。

はっきり言って２＊２をガウスの消去法で解く意味が分からんｗ
というわけだから張り切って一般のｎ次行列について作れｗ

43:35
06/01/28 20:29:59
>>41
これをプログラミングで解いてレポートして提出しろという課題なんですよ。
一応実際に計算して解くことは出来ました。
しかしプログラミングのほうが分からず…
>>42
ガウスは使わない方がいいのでしょうか。
となると逆行列を用いての計算のほうがやりやすいのでしょうか？

44:デフォルトの名無しさん
06/01/28 20:32:57
>>43
じゃあさ、面積を求める部分をSimpson法やいろんな数値積分法で組んでしまうのはどう？
面積を求める部分をサブルーチンか関数にしてaを渡せばそれなりにプログラミングの形になるし。

45:42
06/01/28 20:49:16
>>43
ガウスの消去法はアルゴリズムが単純なので、小さめの連立方程式を解くための
練習問題としてはよく使われる。行列がsingularな時や誤差が溜まりやすい状況を、
あまり考慮しなければ、比較的楽に作れるだろう。

素朴には、対角要素が全部０で無いなら、１列目からｎ列目まで、順に各列の対角要素以外を
０にするように行ごとひいていけばいい。手でとくのと同じだ。
対角要素が０や０付近だと例外処理をすることになる。このへんがまんどくさい。

この辺は、問題製作者の意図がどの辺にあるのかで、こだわり度が決まってくるのだが、
この問題設定では、いまいち理解しにくい。

２＊２での一般形なら単純に公式を使うほうが早いが、n*nの一般サブルーチンを求めているなら、
ガウス法なり何らかの方法が必要になる。

後半では数値積分をすることになるわけだが、問題として全く意図が読めないので、
これだけの情報では、どの方法をとるべきか答えがたい。

46:35
06/01/28 21:01:23
>>44
うーん、すいません、またｔくの初心者ですのでよくわかりません…
助言は非常に嬉しいのですが…すいません。
>>45
今回はどうやら2*2行列のようです。
それでは今回は公式の方を用いて計算することにします。
まぁ、普通の行列のプログラミングもよくわからないのですが…
後半部分はあれだけです。
タダ単に前半で求めた値を代入して面積を求めるだけでして、背景は無いと思います。

47:デフォルトの名無しさん
06/01/28 22:46:30
>>35
まあ、適当にやれｗ

a = 1.0, b = 0.0

交点
x = 0.0, x = e - 1

積分　ｆ（ｘ） = ax - x log(a + x) (a = 1.0　積分範囲　[0.0, e-1])

解析的な値
-e^2 / 4 + e - 5 / 4

台形公式　短精度　分割数 n=1000 で　こんな感じだ。　
誤差は( (e-1)/1000 )^2～1e-6位なので、まぁもっともな気がする。

numerical integral 0.378981858
analytic value　　　0.378982186
Program Completed
Press Enter to Continue.

48:デフォルトの名無しさん
06/01/28 22:51:17
PROGRAM vip
IMPLICIT NONE
INTEGER, PARAMETER :: n = 1000
REAL :: x0, x1, s
x0 = 0.0
x1 = EXP(1.0) - 1.0
s = trapez_integ(n, x0, x1)
PRINT *, 'numerical integral', s
PRINT *, 'analytical value ', 5.0 / 4.0 + EXP(1.0)**2 / 4.0 - EXP(1.0)
STOP
CONTAINS
REAL FUNCTION trapez_integ(n, x0, x1)
INTEGER, INTENT(IN) :: n
REAL , INTENT(IN) :: x0, x1
INTEGER :: i
REAL :: s, h
h = (x1 - x0) / REAL(n)
s = 0.0
DO i = 1, n
s = s + 0.5 * ( func(x0 + (i - 1) * h) + func(x0 + i * h) )
END DO
trapez_integ = s * h
RETURN
END FUNCTION trapez_integ
!
REAL FUNCTION func(x)
REAL, INTENT(IN) :: x
REAL, PARAMETER :: a = 1.0
func = a * x - x * LOG(a + x)
RETURN
END FUNCTION func
END PROGRAM vip

49:デフォルトの名無しさん
06/01/28 22:56:38
訂正

解析的な値
-e^2 / 4 + e - 5 / 4
↓
e^2 / 4 - e + 5 / 4

正負逆だったｗ

交点が　0.0 と　e-1　なのは　func(0.0) と func(e-1) を書かせて確かめるのもよし。

前半は、２＊２の逆行列の公式でいいだろう。
Ax＝ｂ　
ｘ＝（A^-1）ｂ
で一応行列の積も入るし。

50:デフォルトの名無しさん
06/01/29 12:21:03
ひょっとして、数値計算プログラム作るなら
Cよりフォートランの方が高速コードを生成しますか？

51:デフォルトの名無しさん
06/01/29 13:06:57
>>50
昔はそうだった。最近はそれほどで差が無い。
コンパイラ自体がフロントエンド部とバックエンド部に分けて作られるようになり、
中間言語→機械語生成のバックエンド部は共通になったから。

しかしF90以降ならどう考えてもFORTRANの方がバグが少なく、
簡潔明快なソースでやりたいことを表現できる。

52:デフォルトの名無しさん
06/01/29 13:12:43
速度は？

53:デフォルトの名無しさん
06/01/29 13:26:03
スパコンのレベルになると、事実上FORTRANしか選択肢がなかったりする

でもそこらのパソコンでやるなら、Cで書こうがFORTRANで書こうが大して変わらんね。
g77とかはほとんど内部でf2cを動かして、C言語に変換してからgccでコンパイルしてるし

54:デフォルトの名無しさん
06/01/29 13:44:16
>>52
コンパイラメーカーが同じなら、速度もたいして変わらない。
要するに、ソースを中間言語に直すところでの最適化の違いだけ。
中間言語から機械語へは共通化している。
ただ、FORTRANの方が制約がきつい分コンパイラーの判断できる最適化は
Cより高いとされている。

>>53
スパコンに関して言えば、９０年代はじめごろまではベクトル化コンパイラーはFORTRAN
だけで、Cコンパイラーがあっても非ベクトル型だった。しかし、その後Cコンパイラーも
ベクトル化コードを吐くようになってきた。

55:デフォルトの名無しさん
06/01/29 14:03:46
FORTRAN９０／９５はコードを書いていて気持ちいい。
柔軟性と制約のバランスが取れていてぴったりくる手袋や靴下のよう。

イラつくのは宣言部が重くなりがちな点くらいか。

みなにお勧めする。

56:デフォルトの名無しさん
06/01/29 14:09:26
じゃあ、特に速度を要求する部分を、書けるならアセンブリで、書けないならフォートランで、
んで普段はC++で書くとかって感じですかね

57:デフォルトの名無しさん
06/01/29 14:26:58
MPIのルーチンはCとFortranしか提供されていないが、微妙にCより。

C++は動的バインディングのせいでCにくらべて１０～２０％遅くなるという話だ。
次期FORTRANにOOP機能が付け加わったらどうなることか。

58:デフォルトの名無しさん
06/01/29 19:40:26
IMPLICIT REAL*8 (A-H,O-Z)
WRITE(*,*) 'あ'
N=3
COMPLEX*16 A
END

でコンパイルしようとして

$ g77 test3.f
test3.f: In program `MAIN__':
test3.f:2:
WRITE(*,*) 'あ'
1
test3.f:4: (continued):
COMPLEX*16 A
2
Statement at (2) invalid in context established by statement at (1)

って言われるんですけど、何がいけないんでしょう？
教えてください

59:デフォルトの名無しさん
06/01/29 20:22:25
>>58
COMLEX文は宣言文だから、最初の実行文よりも前に置かれねばならない。
このばあい、WRITE文の前。

それで怒られている。

60:デフォルトの名無しさん
06/01/30 10:41:15
2006年の1月から3月までのカレンダーを作れという課題がでました。

・・・・・・・・・・・
do 20 t=1,3
if(t.eq.2)then
b=28
else
b=31
end if
・・・・・・・・・
・・・・・・・・・
write(6,*)" ",t,"月　　　　　　　　"
write(6,*)"　日　月　火　水　木　金　土　　"
・・・・・・・・・・
・・・・・・・・・・
write(6,100)(a(i),i=6,b)
100 format(7i5)
・・・・・・・・・・
・・・・・・・・・・

・・・・の部分がわかりません。よろしくお願いします。

61:デフォルトの名無しさん
06/01/30 11:08:08
いまどきフォートラン使ってる人って、どういう立場だなんだ。

趣味？　仕事？　学校？　環境は？　作ってるアプリの種類は？

ちと興味ある。

62:デフォルトの名無しさん
06/01/30 11:24:39
>>61
失敬な！
FORTRANは現役ですよ！いまやモジュール型言語として完成。次はOOPｗ
INTELが直接開発・販売しているコンパイラはC＋＋とFORTRANだけですよん。

63:デフォルトの名無しさん
06/01/30 12:37:58
いやだから何作ってるのか教えてよ。現役なのは存じてるから

64:デフォルトの名無しさん
06/01/30 15:13:53
FORTRANを科学技術計算以外で使う人なんているの？

65:デフォルトの名無しさん
06/01/30 15:35:34
ということはここに居る人みんな科学者！？

66:デフォルトの名無しさん
06/01/30 16:43:06
>>61, >>63
俺は固体の量子計算。学校の研究です。
アプリっていうとGUIを想像してしまうけど、そうじゃなくて、
コマンドラインで実行するようなプログラム。

>>65
たぶんそうだと思う。
あとガッコウの宿題ができない人々…。

67:デフォルトの名無しさん
06/01/30 17:29:24
Fortran使って構造計算してる土建屋も科学者か？

68:デフォルトの名無しさん
06/01/30 21:27:22
FORTRAN　ｍｐ３　ｐｌａｙｅｒ
URLﾘﾝｸ(www.pbase.com)

ﾜﾛｽ

69:デフォルトの名無しさん
06/01/30 21:55:33
fortran mp3 encoder
URLﾘﾝｸ(members.at.infoseek.co.jp)

70:デフォルトの名無しさん
06/01/30 21:56:14
ところで、fortranに高階関数はいつ入るんですか？

71:デフォルトの名無しさん
06/01/30 22:25:00
>>70
高階かんすうってなんですか？おしえてかださい

72:デフォルトの名無しさん
06/01/31 01:46:45
　　　　　　　　　　　 ,,.　 -‐―-　､
　　　　　　　　　／ _..-‐―- 、　　.ヽ
　　　　　　　　 / ／-‐―-..__..ヽ　_..、ヽ
　　　　　　　　.!三三三´　　　｀´　　:　i
　　　　　　　　|.!三三　　 __　　￣ _,. .;. i　暗黙の型宣言を守らない
　　　　　　　　|ii　i!"　　´ ｪｪ｀　l´rｪ､ i..!　　　　　　>>60は地獄へ逝く
　　　　　　　　 !ｉ､ i!:::..... -￣r　　ヽ　　.!
　　　　　　　　　!ｉｉ!.|:::::::　　/ ｀－´ヽ　/
　 _,,..　===＝"/%;;;::::.　　　'∠ﾆﾆ」　 /、　　　　　　　　　　　　_　-‐‐-､
/´　　　　　　 | .%　:::;;;;:　,,,.,,,:. ⌒: ./´　￣＝=-- 　_,,.. -‐‐'"´　_,..-‐ ´
　　　　　　　 |　 %　　::::::;;;;;;;;;:.""%.／　　　　　　ｒ"　(<二￣￣>
　　　　　　　.|　　%｡ﾟ+。　　。+% ／　　　　　　　｀ヽiﾉ<二.￣￣>
　　　　　　　|　　　%+*ﾟ+。*ﾟ%／　　　　　　　　　（　<二￣￣>

73:デフォルトの名無しさん
06/01/31 02:26:28
>>60　特別にFORTRAN77で書いてやったです。　感謝しやがれです。
PROGRAM calndr
C calender for Jan to Mar 2006
INTEGER ia(31)
CHARACTER*30 fmt
mon1 = 1
C Jan 1 2006 is Sunday (= 1)
DO 10 i = 1, 31
ia(i) = i
10 CONTINUE
DO 20 it = 1, 3
IF (it .eq. 2) THEN
ib=28
ELSE
ib=31
END IF
WRITE(fmt, '(a, i2, a, i1, a)')
1 '(', 4 * (mon1 - 1), 'x, ', (8 - mon1), 'i4, / (7i4))'
mon1 = MOD(ib + mon1 - 1, 7) + 1
WRITE(6, *) ' ', it, '月　　　　　　　　'
WRITE(6, *) '　日　月　火　水　木　金　土　　'
WRITE(6, fmt) (ia(i), i = 1, ib)
WRITE(6, *)
20 CONTINUE
STOP
END

>>73

74:73
06/01/31 02:28:02
実行結果です。

1 月　　　　　　
　日　月　火　水　木　金　土
1 2 3 4 5 6 7
8 9 10 11 12 13 14
15 16 17 18 19 20 21
22 23 24 25 26 27 28
29 30 31

2 月　　　　　　
　日　月　火　水　木　金　土
1 2 3 4
5 6 7 8 9 10 11
12 13 14 15 16 17 18
19 20 21 22 23 24 25
26 27 28

3 月　　　　　　
　日　月　火　水　木　金　土
1 2 3 4
5 6 7 8 9 10 11
12 13 14 15 16 17 18
19 20 21 22 23 24 25
26 27 28 29 30 31

Press any key to continue

75:73
06/01/31 02:31:04
画面が崩れたので、専ブラ用のアンカーを書いてやるです。　
>>73,75

FORMATの技巧が高度すぎて、>>60はもちろん>>60に出題した
どんくさいチビ人間には理解できないだろです！

76:73
06/01/31 02:31:53
専ブラ用のアンカー間違ったです。　
>>73,74

77:デフォルトの名無しさん
06/01/31 03:14:40
>>71 関数を値として引数に取れる関数だよ。

78:デフォルトの名無しさん
06/01/31 07:00:31
関数のポインタ？

79:デフォルトの名無しさん
06/01/31 09:50:36
>>77
単にexternalなら今でも出来るが、関数ポインターならFORTRAN2003で可能に。

80:デフォルトの名無しさん
06/01/31 16:21:20
今日の宿題はまだか？ｗ

81:デフォルトの名無しさん
06/01/31 22:44:21
ところで型推論っていつFortranにはいるんですか？

82:デフォルトの名無しさん
06/01/31 23:40:14
>>81
型推論ってなんですか？

83:デフォルトの名無しさん
06/02/01 00:18:53
すみません、質問なのですが。
フォートランで作った実行ファイルとソースファイルから
そのソースでcallから呼び出されてる関数
のライブラリ(モジュール？）の場所を探す方法はありますか？

というのも、プログラムを使っていたＵＮＩＸマシンから別のＵＮＩＸマシンにプログラムを移したところ
マルチバイト文字が正しくありません。とエラーが出てしまいます。
なのでソースがあるので別のマシンでコンパイル、リンクしようと思ったのです。
そうするとソースの中にcallで呼び出されたサブルーチン記述されていない名前の関数が
あったのでリンクのときにこのモジュールを指定しなければいけないんですよね？
その場所がわからずコンパイルでエラーが出てしまいます。

長文で申し訳ないのですがよろしくお願いします。

84:デフォルトの名無しさん
06/02/01 00:56:50
>>83
エラーメッセージを見ていないので状況がいまいちつかめないが、
MAP出力とかのコンパイラ（リンカ、ローダー）オプションで、
必要な情報が得られるかもしれない。

85:83
06/02/01 20:00:59
>>84
レスありがとうです。

いろいろ悩みましたが
ｺﾏﾝﾄﾞ　nm を使ってライブラリを全部grepに渡して関数名を検索したら出てきました。
新しいマシンに移してコンパイルリンクでき、実行も出来ました。

86:デフォルトの名無しさん
06/02/01 22:34:43
今日は宿題無いのか？

87:こんにちは
06/02/01 22:44:03
日付を2つ入力して，その間に何日あるかを書き出すプログラム(Nannichi)を作れ。実行例はつぎの通りである．日付を２つ入力して下さい．
何年？ 2004
何月？ 7
何日？ 1
何年？ 2005
何月？ 7
何日？ 1
日数は 365日です．
これの答えを教えてください　お願いします

88:デフォルトの名無しさん
06/02/01 22:53:21
>>87
>>5-6をみろ。すでにより一般的なうるう年を考慮したサブルーチンが存在している。
それをちょっといじればよい。

質問があれば、日本語でおｋｗｗｗｗｗ

89:デフォルトの名無しさん
06/02/02 07:55:46
>>82 型を宣言しなくても勝手に正しく推論してくれるしくみだよ。

いつはいるの？

90:デフォルトの名無しさん
06/02/02 10:02:57
>>89
率直に言ってそれは下らないアイデアだと思う。
書きなぐるスクリプト言語ではよくあるが、なんんというかだらしなさのあらわれだ。

どっちにしろ暗黙の型宣言があるし。

91:デフォルトの名無しさん
06/02/03 01:39:04
ROZEN　MAIDEN ・・・　orz

92:デフォルトの名無しさん
06/02/03 23:41:27
ところで、Fortranにproof-carring-codeはいつ入るの？

93:デフォルトの名無しさん
06/02/03 23:49:39
>>92
proof-carring-codeってなんでちゅか？ばぶー

94:デフォルトの名無しさん
06/02/04 00:46:34
だれかこのプログラミングを教えてください。お願いします。
ｘ=0.1 , 0.2 , 0.3 , 0.4 , 0.5
ｙ=1.228 , 1.005 , 0.823 , 0.674 , 0.552
近似式　ｙ=a*exp(bx)
上の関数とデータで最小自乗法で近似した式を求めよ。
この問題のプログラムがわかりません・・・。お願いします。

95:デフォルトの名無しさん
06/02/04 01:35:39
>>94
PROGRAM vip
IMPLICIT NONE
INTEGER, PARAMETER :: n = 5
REAL, PARAMETER :: x0(n) = (/0.1 , 0.2 , 0.3 , 0.4 , 0.5 /), &
y0(n) = (/1.228, 1.005, 0.823, 0.674, 0.552/)
REAL :: y(n), t(n, 2), z(2, 2), zinv(2, 2), v(2)
INTEGER :: i
y = LOG(y0)
t(:, 1) = 1.0
t(:, 2) = x0
z = MATMUL(TRANSPOSE(t), t)
v = MATMUL(TRANSPOSE(t), y)
zinv = ainverse(z)
v = MATMUL(zinv, v)
WRITE(*, *) ' least square fit a=', EXP(v(1)), ' b=', v(2)
STOP
CONTAINS
FUNCTION ainverse(x)
IMPLICIT NONE
REAL, INTENT(IN) :: x(2, 2)
REAL :: ainverse(2, 2)
REAL :: determinant
determinant = x(1, 1) * x(2, 2) - x(1, 2) * x(2, 1)
IF (ABS(determinant) < 1.0e-7) STOP ' singular matrix '
ainverse(1, 1) = x(2, 2) / determinant
ainverse(2, 2) = x(1, 1) / determinant
ainverse(1, 2) = -x(1, 2) / determinant
ainverse(2, 1) = -x(2, 1) / determinant
RETURN
END FUNCTION ainverse
END PROGRAM vip

96:デフォルトの名無しさん
06/02/04 01:36:42
>>95
least square fit a= 1.499271 b= -1.998701
1.227659 1.005252 0.8231379 0.6740159 0.5519093
1.228000 1.005000 0.8230000 0.6740000 0.5520000

97:94
06/02/04 03:35:31
>>95さん　>>96さん
ほんとありがとうございます。
月曜にでも学校でやってみます。
また、わからないことがあれば質問させてください。
ありがとうございました。

98:デフォルトの名無しさん
06/02/04 18:37:30
>>95について。
y=a*exp(bx)を次の式に変形し
log(y)=log(a)+bx
yとxの各値を上の式に代入すると
次の5本の式となる。
log(y1)=log(a)*1+b*x1
log(y2)=log(a)*1+b*x2
log(y3)=log(a)*1+b*x3
log(y4)=log(a)*1+b*x4
log(y5)=log(a)*1+b*x5
ここで、
log(yi)がi番目の要素である列ベクトルをy、
i行1列を1、i行2列をxiとする行列をt、
1番目の要素をlog(a)、2番目の要素をbとする
列ベクトルをvとすると、上の5式は、
y=tv
と表せる。
tが1次独立な列から成るので、
vの最小二乗推定は、次の式
(t't)**(-1)t'y
となる定理を>>95は使用していると思います。
(注)t'はtの転置を表す。(t't)**(-1)はt'tの逆行列を表す。）

99:デフォルトの名無しさん
06/02/04 18:56:37
ところでＦｏｒｔｒａｎに配列のバウンダリチェックはいつ入るのでしょうか？

100:デフォルトの名無しさん
06/02/04 19:05:46
>>99
subscript check は、debug 用にoptionで入るのが普通。1960年代からある。
しかし、遅くなるから実際の実行時に使う奴はいない。
この辺はスピード優先の数値計算の価値観を反映している。

101:９５
06/02/04 19:40:40
オッス！おれ極右！
>>98氏解説㌧

最近寒いので湯たんぽ代わりにノートパソコンを腹の上において寝ているお。
昨夜そういう寝そべった状態で布団の中で作ったのが>>95だお。
オラはVipperだから間違っていても謝罪も賠償もしないお！

最小二乗法なんて普段使わないから、きょうの昼おっきしてからNumericalRecipesで
式を確かめといたおｗ　方式はちょっと違うが式おｋｗｗｗ　（15.2.4~6)あたり。
URLﾘﾝｸ(www.library.cornell.edu)

FORTRAN77で書くと３倍くらいの長さになるとおもうお。
MATMULやTRANSPOSEはF90にしかない行列積と転置だお。
LOGを取ったのは線形化するためだお。
普通、最小二乗法は線形の式にしか使わないと思ったから、そうしてみたおｗ

NumericalRecipesは、旧版のC言語版なら訳本が出ているお。
URLﾘﾝｸ(www.amazon.co.jp)
本のコードにはバグがちりばめられているので、本を打ち込むのはお勧めしないお。
金出して買うコードのほうはこまごまと修正されているお。

（＾ω＾）このくらいソースコードが短いと勉強する気にもなるだろ？

102:デフォルトの名無しさん
06/02/05 07:53:27
>>87
フォートランはわかんないからベーシック
あってるかどうかもわかんない
こんにちは

10 gosub 60
20 nannichi0=nannichi
30 gosub 60
40 print "日数は";nannichi-nannichi0
50 end
60 input "年,月,日",year,month,day
70 if month<3 then year=year-1:month=month+12
80 nannichi=year+int(year/4)-int(year/100)+int(year/400)+int((153*month-457)/5)+day
90 return

103:デフォルトの名無しさん
06/02/05 11:38:27
>>102
解説頼む

104:デフォルトの名無しさん
06/02/05 20:01:35
>>103
まちがいた。ちょと訂正。

80 nannichi=365*year+int(year/4)-int(year/100)+int(year/400)+int(306*(month+1)/10)+day-428

「これは Fairfield公式といって、西暦1年1月1日からの経過日数を算出する式らしい。
　ただしグレゴリオ暦が採用されたのは1582年10月15日からだが。
　要するにこれを使って2つの日付の差を求めればよい。」
っておばあちゃんがいってた。

105:94
06/02/05 20:12:40
すいませんがこないだのプログラムを
ＤＯ文で作っていただけませんか。
何回もお願いして、すいません。

106:デフォルトの名無しさん
06/02/05 20:33:52
コンピューターおばあちゃん♪
コンピューターおばあちゃん♪
イエーイ　イエーイ　僕は大好きさ～♪

107:デフォルトの名無しさん
06/02/05 21:07:36
>>104
把握したｗ
２月の２８日の例外を２月をお尻にすることで回避。
大の月と小の月の入れ替えを１ヶ月の平均日数を30.6日にすることで実現。
かなり絶妙な公式だ。ｵﾓｼﾛｽｗ

だが最近流行の明示的なプログラミングから離れているがｗｗｗｗｗ

108:デフォルトの名無しさん
06/02/05 21:54:46
>>105　世話が焼けるッス！　FORT７７で、素朴に書いといたからｗ
PROGRAM vipper
　　　PARAMETER (n = 5)
REAL x(n), y(n)
　　　DATA x /0.1 , 0.2 , 0.3 , 0.4 , 0.5 /
DATA y /1.228, 1.005, 0.823, 0.674, 0.552/
　　　sumx = 0.0
　　　sumy = 0.0
　　　sumxx = 0.0
　　　sumxy = 0.0
sum = REAL(n)
DO 10 i = 1, n
　　　 sumx = sumx + x(i)
sumy = sumy + LOG(y(i))
sumxx = sumxx + x(i) * x(i)
sumxy = sumxy + x(i) * LOG(y(i))
10 CONTINUE
det = sum * sumxx - sumx * sumx
alog = ( sumxx * sumy - sumx * sumxy ) / det
b = ( sum * sumxy - sumx * sumy ) / det
WRITE(6, *) ' a = ', EXP(alog), ' b = ', b
　　　END

>>108　インデント

109:やり直しｗ
06/02/05 22:11:27
PROGRAM vipper
PARAMETER (n = 5)
REAL x(n), y(n)
DATA x /0.1 , 0.2 , 0.3 , 0.4 , 0.5 /
DATA y /1.228, 1.005, 0.823, 0.674, 0.552/
sumx = 0.0
sumy = 0.0
sumxx = 0.0
sumxy = 0.0
sum = REAL(n)
DO 10 i = 1, n
sumx = sumx + x(i)
sumy = sumy + LOG(y(i))
sumxx = sumxx + x(i) * x(i)
sumxy = sumxy + x(i) * LOG(y(i))
10 CONTINUE
det = sum * sumxx - sumx * sumx
alog = ( sumxx * sumy - sumx * sumxy ) / det
b = ( sum * sumxy - sumx * sumy ) / det
WRITE(6, *) ' a = ', EXP(alog), ' b = ', b
C check
WRITE(6, *) ' i x y',
1 ' a * EXP(bx) (dy)^2'
DO 20 i = 1, n
yy = EXP(alog + b * x(i))
WRITE(6, *) i, x(i), y(i), yy, (y(i) - yy)**2
20 CONTINUE
STOP
END

>>109　インデント

110:デフォルトの名無しさん
06/02/05 22:12:16
a = 1.499271 b = -1.998701
i x y a * EXP(bx) (dy)^2
1 0.1000000 1.228000 1.227659 1.1664589E-07
2 0.2000000 1.005000 1.005252 6.3688631E-08
3 0.3000000 0.8230000 0.8231379 1.9023346E-08
4 0.4000000 0.6740000 0.6740159 2.5326941E-10
5 0.5000000 0.5520000 0.5519093 8.2298044E-09
Press any key to continue
>>110　インデント

111:デフォルトの名無しさん
06/02/05 22:28:29
>>94
死ぬほど易しく作ったから理解しろｗ
>>101の本の公式通りに作り、変数名も大体あわせておいた。
ここで誤差のσは各要素で共通一定と仮定した。この場合キャンセルしあうから式から無視してよし。

本を見れば分かるが、二乗誤差が極値をとる条件を、a,bに関してそれぞれ微分を取ることで求めている。
二乗誤差は極値を１個だけもちかつそれは最小だからこれでおｋとなる。
（誤差は外せばいくらでも大きく出来るが、合わせるほうはに限界がある。最高に合いまくって０。）

プログラム解説：
１０番LOOPで総和を取る。
その後、２＊２の逆行列をかけるのと同じ理屈で、２元連立方程式を解いている。
２０番LOOPでどのくらい近似できたか見ている。

求めていないが、一番右側の列の二乗差の総和が最小になっているはずである。

a,bの値を手動で変えて、求めた二乗差の総和が最小であることを確かめてみそ。
その結果を考察に書いておけば＋１０点の花丸ハンバーグだ！

112:94
06/02/06 00:15:35
みなさんほんと詳しく教えていただきありがとうございます。
ほんと助かりました。

113:デフォルトの名無しさん
06/02/06 12:36:36
ところで、Fortranにオブジェクトデータ構造はいつ入るんですか．

114:sage
06/02/06 13:32:23
>>113
FORTRAN2003。コンパイラーの実装待ち。

115:デフォルトの名無しさん
06/02/06 17:44:37
>>113
オブジェクトデータ構造ってなんですか？
オブイェークト

116:デフォルトの名無しさん
06/02/06 19:53:17
ERK=ER(K)*DCOS(GANMA2(K))-EI(K)*DSIN(GANMA2(K))
EIK=EI(K)*DCOS(GANMA2(K))+ER(K)*DSIN(GANMA2(K))

実数 ERK=ERK*DCOS(GANMA3(K))-EIK*DSIN(GANMA3(K))
虚数 EIK=EIK*DCOS(GANMA3(K))+ERK*DSIN(GANMA3(K))

ERK=ERK*DCOS(GANMA4(K))-EIK*DSIN(GANMA4(K))
EIK=EIK*DCOS(GANMA4(K))+ERK*DSIN(GANMA4(K))

この式をEXPを用いた複素数型に変換できなくこまってます。
どなたか優しい方おしえていただけませんか（＞＜）

117:デフォルトの名無しさん
06/02/06 21:53:43
>>116
やぁ！やさしいお兄さんだよ！援助してあげるから、今日からパパと呼ぶんだよ！
定義：　
ze = er + i ei
zg = COS(Γ) + i SIN(Γ) = EXP(iΓ)

積：
ze * zg = ( er + i ei ) * ( COS(Γ) + i SIN(Γ) )
　　　　　= ( er * COS(Γ) - ei * SIN(Γ) ) + i ( ei * COS(Γ) + er * SIN(Γ) )

これより
COMPLEX*8 :: ze, zek
ze = (er, ei)
zek = ze * EXP( CMPLX(0.0d0, gamma2) )

これでＯＫさ！

ただ細かいことを言うと、倍精度複素数はＦ７７標準規格に存在しないので、
完全に規格に乗っ取ろうとすると複素変数は使えないのさ！
でも、大概のコンパイラは倍精度複素数があるのでＯＫさ。
ＥＸＰとＣＭＰＬＸは総称名なので、精度は引数に合わせてくれるはず。

でもマニュアルを調べて総称名が、ちゃんと倍精度複素数になってるかどうか
確認したほうがいいかもね。

118:117
06/02/06 22:00:01
おっと
ze = CMPLX(er, ei)
だったね。

さぁ、>>116　おっぱいうｐうｐ
　　　 _ 　∩
　　(　ﾟ∀ﾟ)彡　おっぱい!おっぱい!
　　(　 ⊂彡
　　|　　　|　
　　し ⌒Ｊ

119:デフォルトの名無しさん
06/02/06 22:39:00
なんだかしりませんが
Fortran2003って普通にかな～りすごそうなんですけど
ほんとに実装されるんですか？

120:デフォルトの名無しさん
06/02/06 23:17:41
>>118
つ URLﾘﾝｸ(kahouweb.hp.infoseek.co.jp)

121:デフォルトの名無しさん
06/02/06 23:29:50
>>119
FORTRANの場合スーパーコンピュータやCPUのキラーアプリとして
単独での採算を度外視されて開発されているので、多分実装されるでしょう。

INTELが旧DECのFORTRANコンパイラーチームをCompaqから買い取ったのは
そういう方向からだし。

122:デフォルトの名無しさん
06/02/07 11:50:47
次の宿題もってこい！

123:デフォルトの名無しさん
06/02/07 15:28:52
ｆ９０を使い初めてプログラムの複素数化をしているのですがうまくいきません、どなたか助言をおねがいします。
PROGRAM TTT
IMPLICIT NONE
INTEGER :: I,N,M
REAL(8) :: EI,ER(0:2047),EI(0:2047),A,B
REAL(8) :: FACTR,H,PI,RN1,RN2,EP0,C,WL,TINTY,TH,S(0:2047),TT(0:2047),GEFF,ERI,EII
!
PI=4.0*DATAN(DBLE(1.0))
C=3.0E+8
WL=1.55E-6
DO I=0,N-1
FACTR=-H*PI*RN2*RN1*EP0*C/WL
TINTY=(ABS(S(I)*S(I))+ABS(TT(I)*TT(I)))/DBLE(2.0)
TH=FACTR*TINTY
ERI=EXP(-GEFF*H)*(1*DCOS(TH)-1*DSIN(TH))
EII=EXP(-GEFF*H)*(1*DCOS(TH)+1*DSIN(TH))
ER(I)=ERI
EI(I)=EII
OPEN(20,FILE='TH.dat')
WRITE(20,*)ER(I),EI(I)
END DO
END PROGRAM TTT
上のプログラムからSINやCOSをけしてEXP等で表現
し複素数でのプログラムにするにはどうしたらいいでしょうか？

124:デフォルトの名無しさん
06/02/07 16:48:19
>>123
まずその前に・・・・・
C=3.0E+8
WL=1.55E-6
これだと単制度で定数が入るので、ここで数値が腐ってしまう。E→Dにせよ。

吸収係数かなんかの計算か？
発想が逆さまだ。もともとの公式を素直にプログラムするほうがよっぽど楽だぞｗ

ほいで
(1 + i) * EXP(iθ) = (1 + i) * (cosθ+ isinθ)
　　　　　　　　　　= (cosθ+ isinθ) + (icosθ - sinθ)
　　　　　　　　　　= (cosθ - sinθ) + i(cosθ +sinθ)
ゆえに
ZE = EXP(-geff * h) * (1.0d0, 1.0d0) * EXP(CMPLX(0.0d0, th))
こうじゃねーか?

間違っていてもシラネｗ

125:デフォルトの名無しさん
06/02/07 17:30:02
>>124
サンクスです、再度考えてみます

126:124
06/02/07 17:41:11
補足
(1+i) = SQRT(i) * SQRT(2) = EXP( iπ / 2 )

(i + i) * EXP(iθ) = EXP( i(θ + π/2) )

127:124
06/02/07 17:44:09
ごめｗ　今借り物でFEPが違くて勝手が違うｗ　間違い膜リング
(1+i) = SQRT(i) * SQRT(2) = EXP( iπ / 2 )　* SQRT(2)
(1 + i) * EXP(iθ) = SQRT(2) * EXP( i(θ + π/2) )

つーかやっぱ、一度元の公式に戻れｗ

128:124
06/02/07 17:45:34
糞ごめｗ　π/2 →　π/4 うっはーｗｗｗｗｗ　信用度０ｗｗｗｗ

129:124
06/02/07 17:52:55
ZE = EXP(-geff * h) * (1.0d0, 1.0d0) * EXP(CMPLX(0.0d0, th))
　 = EXP(CMPLX(-geff * h, th + pi / 4.0d0)
こうかな？ｗ

スレ汚し　すまんこﾉｼ！

130:123
06/02/07 17:57:00
>>124
できましたありがとうございますm(__)m

131:デフォルトの名無しさん
06/02/07 18:00:55
さぁ、>>123　おっぱいうｐうｐ
　　　 _ 　∩
　　(　ﾟ∀ﾟ)彡　おっぱい!おっぱい!
　　(　 ⊂彡
　　|　　　|　
　　し ⌒Ｊ

132:デフォルトの名無しさん
06/02/07 18:11:54
>>131
つ URLﾘﾝｸ(members.at.infoseek.co.jp)

133:123
06/02/07 18:26:33
実数部分 ERK=ER(K)*DCOS(GANMA2(K))-EI(K)*DSIN(GANMA2(K))
虚数部分 EIK=EI(K)*DCOS(GANMA2(K))+ER(K)*DSIN(GANMA2(K))

ERK=ERK*DCOS(GANMA3(K))-EIK*DSIN(GANMA3(K))
EIK=EIK*DCOS(GANMA3(K))+ERK*DSIN(GANMA3(K))

ERK=ERK*DCOS(GANMA4(K))-EIK*DSIN(GANMA4(K))
EIK=EIK*DCOS(GANMA4(K))+ERK*DSIN(GANMA4(K))
を同じ考えで
EK=E(K)*EXP(CMPLX(0.0d0,(GANMA2(K)+GANMA3(K)+GANMA4(K))))
としたら同結果にならないのですがご教授おねがいします
本当に連投すみませんm(__)m

134:デフォルトの名無しさん
06/02/07 23:14:45
>>133
EK = E(K) * (1.0d0, 1.0d0)**3 * EXP( CMPLX(0.0d0, GAMMA2 + GAMMA3 + GAMMA4) )
= E(K) * SQRT(8.0d0) * EXP(CMPLX(0.0d0, GAMMA2 + GAMMA3 + GAMMA4 + 3.0d0 * pi / 4.0d0) )

でね？

135:デフォルトの名無しさん
06/02/08 01:32:58
>>133はもう一度>>124を読み直しなさい！！ｗ

136:デフォルトの名無しさん
06/02/08 13:12:44
誰かネタを投下しろです！

137:123
06/02/08 16:59:13
前半の距離h/2にわたる線形空間伝搬
[exp( j (r2 + r3 + r4 ) h/2 }F{Φ(z,t)}]
　
元の式がこれなのですが複素数でやろうとするとうまくいかないです・・・
(1.0d0, 1.0d0)の部分に変数をいれればいいのですかね？

138:デフォルトの名無しさん
06/02/08 17:07:52
>>137
その式のままプログラムすればいい。
phi = phi_func(z, t)
zf = F(phi)
zres = EXP( CMPLX（0.0d0, 0.5d0 * h * (r2 + r3 + r4) ) ) * zf

3行に分けたのはdebug時に、１行ごとにエラーを知るため。
うまく行くのがはっきりしたら、１行にまとめればいい。

139:デフォルトの名無しさん
06/02/08 18:12:41
>>138
コンパイラによってはｆｕｎｃｔｉｏｎを別のfunctionやサブルーチンの引数に取れないかもな。

140:デフォルトの名無しさん
06/02/08 22:38:28
1 名前：ぶつわよ！φ ★[] 投稿日：2006/02/08(水) 18:26:54 ID:???0
毎日新聞「開かれた新聞」委員会の１２、１月度見解を報告します。
今回は、昨年１２月２５日に山形県庄内町で起きたＪＲ羽越線の特急脱線転覆事故を扱った
毎日新聞の社説「安全管理で浮ついてないか」（１２月２７日朝刊）が
「風の息づかいを感じていれば、事前に気配を感じていたはずだ」と指摘したことなどに対し、
読者から多数の批判意見が寄せられた問題を中心に取り上げます。

※とっても長いので全文はソースで
URLﾘﾝｸ(www.mainichi-msn.co.jp)

141:デフォルトの名無しさん
06/02/09 01:12:17
プログラミングの宿題です。助けて下さい。

2種の生物x(t)とy(t)が同一場所に生息している。
x(t)の出生率はa>0であり、x(t)はy(t)に食べられる、その死亡率はby(t)<0であるとする。
一方x(t)を食べるy(t)の死亡率はc<0であり、食料となるx(t)に比例して出生率がdx(t)>0であるとする。
このとき次の連立微分方程式系のシミュレーションのプログラムを作れ。
dx/dt=(a+by)x
dy/dt=(c+dx)y
ただし、t=0、x=3、y=1、h=0.2、a=1、b=-1、c=-1、d=1(hは時間tの刻み幅)
　
Fortran77の形式でお願いします。よろしくお願いします。

142:デフォルトの名無しさん
06/02/09 02:18:09
>>141 ほれ！ティンティンうｐねｗ
　微分方程式の積分を最高簡単にしておいた。誤差が溜まって結果怪しい。
　ルンゲ＝クッタとか真面目なのに自分で治しとけｗ

PROGRAM oppai
PARAMETER (h = 0.2, a = 1.0, b = -1.0, c = -1.0, d = 1.0)
t = 0.0
x = 3.0
y = 1.0
DO 10 i = 1, 100
t = i * h
dx = (a + b * y) * x * h
dy = (c + d * x) * y * h
x = MAX(0.0, x + dx)
y = MAX(0.0, y + dy)
WRITE(6, *) x, y, t, i
10 CONTINUE
END

>>142　インデント

結果はというと、ｔ＝１０くらいで、生命体ｘが滅亡し、その後餌のなくなった
生命体ｙもすぐに滅亡したｗ

143:142
06/02/09 02:27:24
>>141
メッシュを細かくして計算したら、ｘとｙは交互に増えたり減ったりして振動したｗ
この結果がもっともなので、>>142のままではやはりﾔｳﾞｧｽｗｗｗｗｗｗｗｗｗｗ

144:コンパイル？
06/02/09 05:24:58
fortran(パソコン？）に初めて触ると言うくらいの初心者です。
URLﾘﾝｸ(www.coastal-env.k.u-tokyo.ac.jp)にそって、".Net"が無い為、fortran77をダウンロードしました。
書いてあるとおり、"F77.exe"の位置を指定して、プログラムをコンパイル(実行）を行った結果、画面には…

「コンパイルできません。
ファイル" …(ファイル名).exe"は、存在しません。」

と表示され、何も反応してくれません。
どうしたらよいでしょうか？教えてください…。

145:コンパイル？
06/02/09 05:50:14
ちなみに、f77.exeは、"C:\Program Files\Fortran77\ftn77.exe"に保存してあります。
保存先は、"…\My Documents\Fortran"に指定しました。

F95でも、同様な事をすれば"fcpad231.lzh"でうまく行くと書いてありますが、こちらも教えて頂けましたら、幸い名のですが…。
(今、自分は研究者なのですが…遅まきながらの数値計算の導入で、かなり四苦八苦しています。よろしくお願いします。）

146:デフォルトの名無しさん
06/02/09 07:03:53
>>144
それは
「あなたの作ったプログラムに間違いがあるのでコンパイルできませんでした。
（コンパイルできなかったので）実行ファイルが存在せず、実行できません」
という意味です。
「メッセージ」の欄にエラーの起こった行とその理由が表示されるので
それを参考にもう一度よくプログラムを見直して訂正して下さい。

>>145
FTN95 コンパイル時のオプション指定の方法が FTN77 と違うため、
そのままでは CPad for Salford FTN77 では使えません。
FTN95付属のIDEを使うか（但し、日本語は使えません）、
適当なエディタを使って下さい。

147:コンパイル？
06/02/09 07:52:35
>>146
ありがとうございます。そうだったんですね。

メッセージには、
■C:\Documents and Settings…\デスクトップ> ftn77 /link ex.for
[CPad: コンパイル終了]
とあります。

どうしたらよいのでしょうか…？

あと、適当なエディタを用いるとありますが、どのようなものを使用すればよいのでしょうか？
本当に初歩的な質問ですいません…。

148:デフォルトの名無しさん
06/02/09 08:03:53
もう一つ質問させてください。
f77よりf90でパンコイルしたほうが高速ですか？

149:デフォルトの名無しさん
06/02/09 11:33:46
>>148
初心者レベルで問題になるような著しい差はない。

エディターに関しては、初心者だというなら

何がしたいのか分からないが、新たに言語を習得したいというならF90を
学ぶほうがいい。

150:デフォルトの名無しさん
06/02/09 11:38:41
エディターに関しては、初心者だというなら予約語（？）が色づけされる
タイプのものがいいはずなので、統合環境のものを使うのが吉だろう。

MSの開発環境は複雑化しすぎて初心者には荷が重いだろう。

無料のものといえば、旧ソ連の金のない連中が、g77や著作権切れの
本のPDFなどをセットにしたFORTRAN統合環境を出していたがｗｗｗ
ま、FTNでいいんでないの？

151:デフォルトの名無しさん
06/02/09 20:29:47
>>147
>どうしたらよいのでしょうか…？

あなたが使っているPCの環境やプログラムの内容がわかるエスパーが現れるまで
しばらくお待ち下さい。

152:デフォルトの名無しさん
06/02/10 02:28:09
>>141　今日は2次のルンゲ=クッタにしてみたｗ
PROGRAM oppai2
PARAMETER (h = 0.2, a = 1.0, b = -1.0, c = -1.0, d = 1.0)
t = 0.0
x = 3.0
y = 1.0
hh = h / 2.0
DO 10 i = 1, 200
t = i * h
dx0 = (a + b * y) * x * hh
dy0 = (c + d * x) * y * hh
dx = (a + b * (y + dy0) ) * (x + dx0) * hh
dy = (c + d * (x + dx0) ) * (y + dy0) * hh
x = MAX(0.0, x + dx)
y = MAX(0.0, y + dy)
WRITE(6, *) x, y, t, i
10 CONTINUE
END

>>152　インデント

とりあえず、振動が続くようになった。捕食のモデルから予想されるものに
なっていると思う。
ま、普通は４次のルンゲ=クッタなんだがｗ

153:デフォルトの名無しさん
06/02/10 02:47:02
>>152　4次ルンゲ食った
PROGRAM oppai4
PARAMETER (h = 0.2, a = 1.0, b = -1.0, c = -1.0, d = 1.0)
t = 0.0
x = 3.0
y = 1.0
hh = h / 2.0
DO 10 i = 1, 2000
t = i * h
dx1 = (a + b * (y ) ) * (x ) * hh
dy1 = (c + d * (x ) ) * (y ) * hh
dx2 = (a + b * (y + dy1) ) * (x + dx1) * hh
dy2 = (c + d * (x + dx1) ) * (y + dy1) * hh
dx3 = (a + b * (y + dy2) ) * (x + dx2) * hh
dy3 = (c + d * (x + dx2) ) * (y + dy2) * hh
dx4 = (a + b * (y + dy3) ) * (x + dx3) * hh
dy4 = (c + d * (x + dx3) ) * (y + dy3) * hh
x = MAX(0.0, x + dx1 / 6.0 + dx2 / 3.0 + dx3 / 3.0 + dx4 / 6.0)
y = MAX(0.0, y + dy1 / 6.0 + dy2 / 3.0 + dy3 / 3.0 + dy4 / 6.0)
WRITE(9, *) x, y, t, i
10 CONTINUE
END

>>153　インデント
振動しつつ、ｘ、ｙとも１に収束してゆく。　ｔ＝３５０　くらいでほぼ収束か。

154:デフォルトの名無しさん
06/02/10 22:57:53
ポトペタ可能な開発環境はありませんか？

155:デフォルトの名無しさん
06/02/10 23:17:14
ポトペタ　ってなんですか？

156:デフォルトの名無しさん
06/02/10 23:34:45
>>154
こんなんがあるな。使ったこと無いけどｗ
URLﾘﾝｸ(www.hulinks.co.jp)

つーか、GUI は他の言語で作った方が楽だけどな。

157:デフォルトの名無しさん
06/02/11 14:41:04
>>147
・・・いっぽうロシア人は鉛筆を使った。・・・・

金なしロシアのｆ７７開発環境だｗ
URLﾘﾝｸ(www.imamod.ru)

158:デフォルトの名無しさん
06/02/11 22:26:04
多変数多項式のGCD計算がなかなか作れないのですが、どなたか見本みたいな感じでソースを見せて頂けませんか？
あと、グレブナ基底を求めるプログラムもあれば頂きたいです。

159:デフォルトの名無しさん
06/02/11 22:48:29
>>158
　【FORTRAN７７】　ﾌﾟﾛｸﾞﾗﾑを起動させて！！
ｽﾚﾘﾝｸ(sim板:68番)

はい、次の方どうぞー

160:デフォルトの名無しさん
06/02/12 01:44:35
>>158　いろいろあった気がするが。　　岩波から本が出てなかったか？　買え！

To search the Netlib repository, enter your search query in the form below and
then press the Go button.

Look for : any of the wordsall the words

Query "groebner" found 1 matches.
1. toms/628
keywords: groebner basis, polynomial ideals, rational integers
gams: C3b
title: GROEB
for: canonical (or Groebner) bases of polynomial ideals
by: F. Winkler et al.
ref: ACM TOMS 11 (1985) 66-78
Score: 100%

URLﾘﾝｸ(www.netlib.org)

161:デフォルトの名無しさん
06/02/12 02:11:34
すごい資料が出たらしいぞｗｗｗｗ

URLﾘﾝｸ(bbs.enjoykorea.jp)

これをみれ！　そして１００回saveしろ！

162:158
06/02/12 02:22:29
>>159
あ～、すみません。それ俺です…
向こうは人が少ないし、偶々ここを見つけて、こっちの方が本職かな？と思いまして。

>>160
岩波からですね。漁ってみます。
自分は脳ミソの出来がよろしくないので、わからないときはまた質問させて頂きます m(_ _)m

163:158
06/02/12 02:39:54
連続ですみません。

>>160
この検索すごいですね。grobnerもGCDも一発で見つかりました。
まだ中身を細かく見てないですが、
岩波の本を買って参照しながら解読してみます。
助かりました。

164:デフォルトの名無しさん
06/02/12 02:55:11
>>162
岩波でなかったかもｗ　間違ってたらすまんｗ
プログラムはよくネットに落ちているので検索すればもっと出てくると思う。

165:デフォルトの名無しさん
06/02/12 02:59:24
>>163
NetLibはインターネット以前からある数学系のサブルーチン集だ。
古いのがごみのように詰まっているｗ　

166:158
06/02/12 04:18:13
>>164
それらしい本で
岩波から計算代数と計算幾何(佐々木,今井,他,共著)
共立から計算機数学(野崎著)
を本屋の検索で発見したので、明日図書館で調べてみようかなと思います
groebnerのプログラムが長すぎて脳がスーパーノヴァしてるので、今日は寝ることにします(_ _;)
数値解析の本ばかりで、GCD計算は代数系の入門書からだったので、本当に助かりました

>>165
古くても重宝しますね、これは。
数にビビリました。
ルンゲも、ガウスも、DE公式も、SVDも、なんでもありますね！
無償で頂けるなんて親切な人達がいるもんですね(^-^)

167:デフォルトの名無しさん
06/02/12 17:15:33
>>166
昔の大型計算センターはどこでもライブラリが大量に公開されているのが普通だった。

ＦＯＲＴＲＡＮの世界は歴史が古いので、ソフト囲い込み運動と
その反動のオープンソース運動とは無関係にソース公開共有の
伝統がまだ残っている。田舎で鍵をかけないで出かける世界？ｗ

168:デフォルトの名無しさん
06/02/14 03:02:05
さぁ次の宿題もってこい！

169:デフォルトの名無しさん
06/02/14 20:32:01
配列を返す関数はどのように定義するのですか？

170:デフォルトの名無しさん
06/02/14 22:20:04
>>169
program hoge
　implicit none
　real, dimension(3) :: x, y, z
　real, dimension(5) :: u, v, w
　x = (/1, 2, 3/); y = (/2, 3, 4/)
　u = (/4, 5, 6, 7, 8/); v = (/10, 20, 30, 40, 50/)
　z = adda(x, y)
　w = adda(u, v)
　print *, z
　print *, w
　stop
contains
　function adda(a, b)
　　real, intent(in), dimension(:) :: a, b
　　real, dimension(size(a)) :: adda
　　integer :: na, nb
　　na = size(a); nb = size(b)
　　if (na > nb) then
　　　adda(1:nb) = a(1:nb) + b(1:nb); adda(nb+1:na) = 0.0
　　else if (na < nb) then
　　　adda = a(1:na) + b(1:na)
　　else
　　　adda = a + b
　　end if
　end function adda
end program hoge

次ページ