[表示 : 全て 最新50 1-99 101- 201- 301- 2chのread.cgiへ]
Update time : 06/21 19:50 / Filesize : 120 KB / Number-of Response : 304
[このスレッドの書き込みを削除する]
[+板 最近立ったスレ&熱いスレ一覧 : +板 最近立ったスレ/記者別一覧] [類似スレッド一覧]


↑キャッシュ検索、類似スレ動作を修正しました、ご迷惑をお掛けしました

くだすれFORTRAN(超初心者用)その2



1 名前:デフォルトの名無しさん [2006/11/22(水) 00:00:36 ]
このスレッドは、他のスレッドでは書き込めない超低レベル、
もしくは質問者自身何が何だか分からない質問を勇気を持って書き込むスレッドです。
FORTRAN使いが優しくコメントを返しますが、
お礼はFORTRANの布教と初心者の救済と次期Fortran2008規格でのCOMEFROM文採用をお願いします。

●注意事項
・質問する前にGoogle等の検索サイトで検索しましょう。
・回答者にわかりやすい様に、質問内容はできる限り詳しく書きましょう。
・エラーの場合は起きた状況、環境(OS・コンパイラ)、エラーメッセージも詳しく書きましょう。


●前スレ
くだすれFORTRAN(超初心者用)
pc8.2ch.net/test/read.cgi/tech/1138063703/

●関連スレ
FORTRAN W
pc8.2ch.net/test/read.cgi/tech/1163319215/


62 名前:デフォルトの名無しさん mailto:sage [2006/12/03(日) 15:22:42 ]
>>47 世話が焼けるなー
PROGRAM vip77
REAL rotz(3, 3), roty(3, 3), vec1(3), vec2(3), tmp(3)
pi = 4.0 * ATAN(1.0)
degrad = pi / 180.0
alpha = 6242.2 / 6370.0
anu = ( 180.0 - 86.6 ) * degrad
theta1 = ( 90.0 - ( 35.0 + 40.0 / 60.0 ) ) * degrad
phi1 = ( 139.0 + 12.0 / 60.0 ) * degrad
C
rotz(1, 1) = COS(phi1)
rotz(1, 2) = -SIN(phi1)
rotz(1, 3) = 0.0
rotz(2, 1) = SIN(phi1)
rotz(2, 2) = COS(phi1)
rotz(2, 3) = 0.0
rotz(3, 1) = 0.0
rotz(3, 2) = 0.0
rotz(3, 3) = 1.0
C
roty(1, 1) = COS(theta1)
roty(1, 2) = 0.0
roty(1, 3) = SIN(theta1)
roty(2, 1) = 0.0
roty(2, 2) = 1.0
roty(2, 3) = 0.0
roty(3, 1) = -SIN(theta1)
roty(3, 2) = 0.0
roty(3, 3) = COS(theta1)



63 名前:デフォルトの名無しさん [2006/12/03(日) 15:23:37 ]
C
vec1(1) = SIN(alpha) * COS(anu)
vec1(2) = SIN(alpha) * SIN(anu)
vec1(3) = COS(alpha)
C
CALL mulmat(3, roty, vec1, tmp)
CALL mulmat(3, rotz, tmp, vec2)
theta2 = ASIN(vec2(3)) / degrad
phi2 = ATAN2(vec2(2), vec2(1)) / degrad
PRINT *, theta2, phi2
STOP
END
!
SUBROUTINE mulmat(n, a, b, c)
REAL a(n, n), b(n), c(n)
DO 10 i = 1, n
c(i) = 0.0
DO 20 j = 1, n
c(i) = c(i) + a(i, j) * b(j)
20 CONTINUE
10 CONTINUE
RETURN
END

64 名前:62 [2006/12/03(日) 16:12:00 ]
使っているのは高校でも習う座標の回転だけだ。
高校では2次元、ここでは3次元だが、y軸とz軸に関して廻すだけだから
高校のものと本質的には変わりはない。COS、SINが湧いて来るのはそれが理由だ。

あと、FORTRANの三角関数はラジアンで与えることになっているから、
度表示をラジアンに直さなければならない。π/180はそこから湧いてくる。

また、角度表示は古代メソポタミア文明以来の60進法が使われているので、
度、分、秒を度の小数点表示にしてから、さらにラジアンにしている。
それが60で割ったりしている部分だ。

あと、座標としては極座標をとるが、北緯、東経などを極座標に直す必要がある。
座標系は右手系を取るのがふつう。(右手系とは右手のお父さん指がx軸、
お母さん指がy軸、のっぽのお兄さん指がz軸だ。)


考え方はこうだ。まず八王子が北極にあるとする。北極からミサイルをぶっ放した
到着点をしるのはたやすい。その地点をvec1に入れている。
ここで北極点を本来の八王子の位置に持っていくように座標軸を回転させる。
このときvec1にも同じ回転をかければ求める地点が分かる。それがvec2だ。

その他、問題文にもあるが、弧の長さと半径が分かれば角度が出る。
方位角ν(読み方はニュー速のニューだぞ)は、本来の八王子から見て真北から
時計回りにとると考える。(問題文からは判然としないが、常識的にはこう取るのが普通)
北極点から考える場合は南極向きで時計回りになる。




65 名前:62 [2006/12/03(日) 16:14:05 ]
訂正

>あと、座標としては極座標をとるが、北緯、東経などを極座標に直す必要がある。

角度の取り方は極座標の取り方をとって、x、y、zの直交座標系を右手系でとった。


66 名前:デフォルトの名無しさん mailto:sage [2006/12/03(日) 20:00:45 ]
>>62
解説乙。
本当は出題者が地球半径を与えている筈ですね。
もし地球を楕円体として考えた場合、かなり面倒そうですね。

67 名前:デフォルトの名無しさん mailto:sage [2006/12/04(月) 02:54:16 ]
>>66
それは強烈に難しいw

68 名前:39 mailto:sage [2006/12/04(月) 10:07:52 ]
>>51
ありがとうございます。出来ました。
subroutineをmoduleで包まないといけないのですね。
最近77から移ってきたのでどうもmoduleの使い方がいまいち
分かりません。

69 名前:デフォルトの名無しさん mailto:sage [2006/12/04(月) 13:54:10 ]
subroutineのなかにsubroutineとかfunction定義できるの?
自分で実験したほうが早いって?www

70 名前:デフォルトの名無しさん [2006/12/04(月) 23:50:42 ]
>>69
できる・・・・が、それは普通のサブルーチンやファンクションとは違っている。
親ルーチンの変数をプチ・グローバルに共有している。

うっかりiとかjとかローカルに宣言しわすれたまま使ってしまうと、
親ルーチンのほうのiやjが書き換わってとんでもないことになる。

好みの問題だが、漏れはあまり好きくない。





71 名前:デフォルトの名無しさん [2006/12/05(火) 14:31:45 ]
CPad for Salford FTN77
のコンパイラのパスが分かりません。
教えて下さい


72 名前:デフォルトの名無しさん [2006/12/05(火) 14:43:44 ]
>>71
>>57のへん見れ

73 名前:デフォルトの名無しさん [2006/12/05(火) 14:51:02 ]
>>72
C:\Program Files\Salford Software\FTN95\ftn95.exe
と入力しても
”コンパイラのパスが正しくありません”と表示されます
ちなみにOSはXPです

74 名前:デフォルトの名無しさん mailto:sage [2006/12/05(火) 14:55:10 ]
>>73
エクスプローラで"ftn95.exe"を検索したら?

75 名前:デフォルトの名無しさん mailto:sage [2006/12/05(火) 22:57:15 ]
>>73
つーかお前の使いたいコンパイラはFTN77かFTN95かどっちなんだ?
FTN77しか入れてないのにftn95.exeを探しても見つかりっこないぞw

76 名前:デフォルトの名無しさん [2006/12/06(水) 15:06:43 ]
FTN77です

77 名前:デフォルトの名無しさん mailto:sage [2006/12/06(水) 15:10:42 ]
もぅ、ぉにぃちゃんしっかりしてよ〜

78 名前:デフォルトの名無しさん [2006/12/06(水) 15:29:18 ]
>>71
なんとかできました

79 名前:デフォルトの名無しさん mailto:sage [2006/12/06(水) 15:42:22 ]
>>78
どうやったら「なんとかでき」たの?

80 名前:デフォルトの名無しさん [2006/12/07(木) 13:24:39 ]
FTN77peをダウンロードして,コンパイラのパスに
FTN77pe〜exeと入力してら、できました





81 名前:デフォルトの名無しさん mailto:sage [2006/12/10(日) 00:22:42 ]
FTN95ってどうよ?


82 名前:デフォルトの名無しさん mailto:sage [2006/12/10(日) 02:43:06 ]
>>81
コンパイル時/実行時のエラーチェックがけっこう強力で
デフォルトで未使用変数の警告なんかも出してくれるので、
入門用やエラーチェック用にイイ感じ。
Personal版がフリーで使えるし。

その他の特徴としては
・Visual Studio 2003/2005上で使える(もちろんVSは別途購入する必要がある)
・VSがなくても独自のIDEが付属する(日本語は使えないが)
・(今時のFortranに必要かどうかはともかく)インラインアセンブラが使える
・(これも実用性のほどはともかく).NETなプログラムも作れる
・GET_COMMAND_ARGUMENT()、[...] 等Fortran2003の一部に対応(他のコンパイラでも大抵使えるレベルだが)



83 名前:デフォルトの名無しさん [2006/12/10(日) 17:07:02 ]
>>82
dクス。
意外に高機能だな。
サポートとかはどうなんだろう?パッチとか結構でてるんだろうか?


84 名前:デフォルトの名無しさん mailto:sage [2006/12/10(日) 20:23:17 ]
intel fortran ver 9.0 for winを使っております。
自分でプログラムを組んでコンパイルできるところまでは到達したのですか実行すると
forrtl: severe (168):Program Exception - illegal instruction
Image PC Routine Line Source

Stack trace terminated adnormally.

と実行してくれません。なにがおかしいのでしょうか?
ソースをのせるべきなんでしょうがnetnumpacから主に引っ張ってきている分権利がややこしそうです。

よろしくお願いします。

85 名前:84 mailto:sage [2006/12/10(日) 20:35:24 ]
program main
common/ numbb / nnn,lll,dz,ameson
integer :: nnn, lll
real*8 :: dz, ameson, dr(251),scr,hx
open (unit=10, file='functest0.txt', status='unknown')
open (unit=11, file='functest1.txt', status='unknown')
open (unit=12, file='functest2.txt', status='unknown')
open (unit=13, file='functest3.txt', status='unknown')
open (unit=14, file='functest4.txt', status='unknown')
open (unit=15, file='functest5.txt', status='unknown')
nnn = 6
lll = 5
dz = 60.0d00
ameson = 413.536d00
dr(1)=exp(-8.8)*dz
hx=5.d-2


86 名前:デフォルトの名無しさん [2006/12/10(日) 20:37:23 ]
do l=2, 251
dr(l)=dr(1)* exp(hx*(l-1))
end do
zero = 0.0
y = 0.0
do i = 1, 251
call GASND (zero, dr(i), screen, 40, Y, ICON)
scr = y
write(10,'("i="i2," dr(i)="e15.8," scr="e20.10)')i,dr(i),scr
end do
close (unit = 10)
close (unit = 11)
close (unit = 12)
close (unit = 13)
close (unit = 14)
close (unit = 15)
end


87 名前:84 [2006/12/10(日) 20:38:49 ]
function screen(ddr)
! it's the Dnl(r)
common/ numbb / nnn,lll,dz,ameson
real*8 :: a1, a2, a3, al1, arn
integer :: nnn, lll, l5
real, external :: DPLAGG
integer, external :: factorial
real, intent(in) :: ddr
a0 = ameson
arn = - a0 * ddr
aa1 = 2.0d00*dz/(nnn*a0)
nl = nnn-lll - 1
nll = nnn + lll
a2 = factorial(nl)
a3 = factorial(nll)
a4 = 2.0d00*dble(nnn)*a3*a3*a3
a5 = aa1*(a2/a4)
l5 = 2*lll + 1
l6 = 2*lll
a6 = (arn)**(l6)
al6 = dble(l5)
al7 = DPLAGG(nll,al6, arn)
al8 = al7*al7
ss = a5*exp(arn)*al8*a6
screen = ss
write(11, '("nnn="i2," ss="e20.10)') nnn, ss
end function screen

88 名前:84 mailto:sage [2006/12/10(日) 20:40:33 ]
function factorial(n)
integer :: fact
integer, intent(in):: n
fact = 1
do i = 1, n ,1
fact = fact*n
end do
factorial = fact
write (12, '("n="i1," fact="i10)') n, fact
end function factorial

FUNCTION DPLAGG (N, A, X)
IMPLICIT REAL (8)(A - H, O - Z)

END FUNCTION DPLAGG
SUBROUTINE GASND (A, B, FUN, N, Y, ICON)
IMPLICIT REAL (8)(A - H, O - Z)
EXTERNAL FUN
DIMENSION X (50), W (50)
SUBROUTINE WGLEGD (NP, PT, WT, ICON)
! LEGENDRE-GAUSS FORMULA
IMPLICIT REAL (8)(A - H, O - Z)
DIMENSION X (25, 50), W (25, 50), PT (NP), WT (NP)


89 名前:84 mailto:sage [2006/12/10(日) 20:42:38 ]
netnumpacからのソースはプログラム名だけにしました
変数の定義があったりなかったりで非常にみづらいとは思いますがよろしくお願いします。

90 名前:デフォルトの名無しさん mailto:sage [2006/12/10(日) 21:08:10 ]
>>84 これは試したかい?
Re: forrtl: severe (168) Program Exception - Illegal instructiom
Reply Quote
Turns out that the problem was caused by an older-generation processor not understanding
newer instructions. The application had been compiled with the "Generate most optimized code"
(/fast) setting, which implies /arch:host. Unfortunately, the project settings box display doesn't
reveal the implied switches, leading to this sort of problem.

Steve

softwareforums.intel.com/ISN/Community/en-US/forums/thread/106764.aspx



91 名前:デフォルトの名無しさん mailto:sage [2006/12/10(日) 21:21:02 ]
>>83
SilverfrostのホームページからForumを見に行くと
ここのAnnouncementsにパッチの情報も載っている。

92 名前:84 mailto:sage [2006/12/10(日) 21:41:42 ]
>>90
英語をよみくだくと /arch:hostのコンパイラオプションをしろ ですよね
やってみましたが変化なしでした。
環境はpen4 2.4BGHz 1024MB 845Echipset winxp sp2なんですがね

93 名前:デフォルトの名無しさん mailto:sage [2006/12/10(日) 22:03:23 ]
>>92
debugモードでtracebackもonにしてどのソース行でおかしくなってるか調べてみたら?


94 名前:デフォルトの名無しさん mailto:sage [2006/12/11(月) 01:23:37 ]
>>92
Pentium 4の2.4BGHzだと拡張命令はSSE2までだから、
どこかにSSE3等を使ったコードが紛れ込んでるんですかね。

別のPCでコンパイルしたライブラリをリンクしているのなら
そのライブラリがSSE3命令を使っている可能性も考えられます。

95 名前:84,92 mailto:sage [2006/12/11(月) 03:06:49 ]
>>93
debugモードは使ったことがないけどやってみます
>>94
ライブラリはもってきていないはずなので
SSE3をつかったコードが紛れ込んでるんですかね


96 名前:デフォルトの名無しさん mailto:sage [2006/12/11(月) 08:32:14 ]
>>95
つーかVisualStudio使わないでコマンドラインでやってるの?デフォはdebugのはずだが。

debugモードを使わないという神経が分からんぞい!
まずは警告オプションも全部onにしてやりんしゃい!

97 名前:デフォルトの名無しさん [2006/12/13(水) 20:37:45 ]
fortran77でファイルに書き込みをするときに
open(2,file='filename.txt')
とかけばfilename.txtに結果が書けるのですが、

ひとつのプログラムでa=1〜100まで変化したときに
(do a=1,100 〜continueを利用)
a=1のときの結果はfilename1.txtに記録
   ・・・
a=100のときの結果はfilename100.txtに記録
するにはどのようにすればいいのでしょうか?

可能であればa=i,b=j(i,jに自然数が入る)のときに
filename_a=i_b=j.txt
に書き込めるように、お願いします。

98 名前:デフォルトの名無しさん mailto:sage [2006/12/13(水) 21:58:54 ]
>>97
前スレでも出てたと思うが、FORTRAN77の場合だと
fmt, fname は適当な長さの文字列変数として

write(fmt,100) int(log10(real(i)))+1, int(log10(real(j)))+1
100 format('(''fname_'',I',I1,',''_'',I',I1,',''.txt'')')
write(fname, fmt) i, j
open(10, file=fname, ...

書式に 'I0' を指定すると「その数値を表現するのに最小の欄幅」を取ってくれる処理系なら

write(fname,'(''filename_'',I0,''_'',I0,''.txt'')') i,j
open(10, file=fname, ...

個人的には数値は適当な桁数にそろえて
filename_001_001.txt, ..., filename_010_010.txt, ..., filename_100_100.txt, ...
としたほうが後でデータファイルをいじったりする時に都合がよいように思う。

write(fname,'(''filename_'',I3.3,''_'',I3.3,''.txt'')') i,j
open(10, file=fname, ...


99 名前:84,92 mailto:sage [2006/12/13(水) 23:06:12 ]
>>96
コンパイラマシンは混んでいていつもコンパイルのみに使っていたのでdebugモードは未経験でした

debugしてみたとこと
dr(1)=exp(-8.8)*dz がひっかかっており
sub = -8.8
dr(1) = dz *exp(sub) としたら動くようになりました。

exp(-8.8)がSSE3なりを使用しているんでしょうか

100 名前:デフォルトの名無しさん mailto:sage [2006/12/14(木) 00:25:28 ]
最適化オプションを無効にすればSSE/2/3は使われないはずだが。



101 名前:デフォルトの名無しさん mailto:sage [2006/12/14(木) 01:14:45 ]
>>99
ReleaseでSSE抑止命令(CPUを旧型に指定)すれば、SSEが原因なのかどうか分かる。

定数式だから変な最適化でバグったのかね?

102 名前:デフォルトの名無しさん [2006/12/14(木) 03:12:02 ]
求めるP3の値がちょうど3(限りなく近い)値になるときのdaの値
を小数点6桁ぐらいまで求めたいのですが、どうすれば次のプログラムを
効率よくdaの値が一回で出せるプログラムにできるか知恵をお貸しください。

DIMENSION E22(2),oo(2)
REAL k,pi,a,o,z,Kf,kf1
pi=3.1415926
DO 30 da= 1.000 , 0.800 , -0.002
DO 30 n4= 8 , 8

k=5.319
a=k/0.529177
o=(a**3)/2.0*da
oo(1)=o
z=1.0


103 名前:デフォルトの名無しさん [2006/12/14(木) 03:12:51 ]
rs=(3.0*o/(z*4.0*pi))**(1.0/3.0)
Kf=(9.0*pi/4.0)**(1.0/3.0)/rs
Ef=Kf**2
c=2.0*1.0/(1.0+0.0155*pi/Kf)
Ei=1.79186*z**(5.0/3.0)/(3.0*o*rs)
P=-Ei*1.4703*10**4
E0=-4.42*z/(3.0*o*rs**2)+0.916/(3.0*o*rs)+0.031/(3.0*o)
P0=-E0*1.4703*10**4
Rmk=3.182
q0=0.89*2*0.3947587
uk=q0*Rmk*cos(q0*Rmk)/(sin(q0*Rmk)-q0*Rmk*cos(q0*Rmk))
E1=-2.0*pi*z*2.0/o**2*(z*Rmk**2*(1.0+2.0/3.0*uk))
P1=-E1*1.4703*10**4
DO 20 L= 1 , 2
E2=0.0

IF(L.EQ.1) GO TO 100
o=o*0.996
oo(L)=o
a=(o*2.0)**(1.0/3.0)
rs=(3.0*o/(z*4.0*pi))**(1.0/3.0)
Kf=(9.0*pi/4.0)**(1.0/3.0)/rs
Ef=Kf**2
c=2.0*1.0/(1.0+(0.0155*pi)/Kf)


104 名前:デフォルトの名無しさん [2006/12/14(木) 03:13:27 ]
100 DO 10 n1= -n4 , n4
DO 10 n2= -n4 , n4
DO 10 n3= -n4 , n4
IF(n1==0.and.n2==0.and.n3==0)GO TO 10
IF(MOD(n1+n2+n3,2)/=0)GO TO 10

G=2.0*pi/a*SQRT(FLOAT(n1*n1+n2*n2+n3*n3))
f=G**2/(2.0*(G**2+c*Kf**2))
ae=1.0+(1.0-f)*4.0*pi*z*2/(o*G**2)*3.0/(2.0*Ef)*
&(1.0/2.0+((4.0*Kf**2-G**2)/(8.0*Kf*G))*ALOG(ABS((2.0*Kf+G)/(2*Kf-G))))
V=(-4.0*pi*Z*2.0/(o*G**2))*((1.0+uk)*cos(G*Rmk)-uk*sin(G*Rmk)/(G*Rmk))

E3=1.0/(4.0*pi*2.0/(o*G**2))*V**2*(1.0/(1.0-f))*(1.0/ae-1.0)
E2 = E2 + E3

10 CONTINUE
E22(L)=E2/2.0
20 CONTINUE

E8=(E22(2)-E22(1))/(oo(2)-oo(1))
P2=-E8*1.4703*10**4
P3=P+P0+P1+P2

WRITE(6,50)P3
50 FORMAT(F15.8)
30 CONTINUE
WRITE(*,*)P,P0,P1,P2
WRITE(6,601) P3
601 FORMAT(1H , F15.8)
STOP
END


105 名前:デフォルトの名無しさん mailto:sage [2006/12/14(木) 10:13:55 ]
>>102
その質問の仕方じゃ 答えようがない。
もう少し何をやっているのか、根本から家

106 名前:97 mailto:sage [2006/12/14(木) 13:14:13 ]
>>98
出来ました。ありがとうございます!

107 名前:デフォルトの名無しさん mailto:sage [2006/12/14(木) 16:39:03 ]
fortran90です。
いまwindows上で動かしてたプログラムを
Linuxに移そうとしてるんですが

Unable to open MODULE file dflib.mod

とエラーが出てきます。
どうすればいいでしょうか?

108 名前:デフォルトの名無しさん mailto:sage [2006/12/14(木) 20:00:26 ]
Unable to open MODULE file dflib.mod
   ↓
エキサイト翻訳
   ↓
MODULEファイルdflib.modを開くことができません。


109 名前:デフォルトの名無しさん [2006/12/14(木) 20:12:39 ]
>>107
ほかの人にやってもらう

110 名前:デフォルトの名無しさん [2006/12/14(木) 22:50:01 ]
>>107
dflib.mod はDECの拡張モジュールだから、そのままではLinuxでは動かない。

DECの拡張ルーチン中で用いているものがPOSIX互換のルーチンだったら
LINUX上に移植できるだろう。

そうでないなら、自分で回避策を考えねばならない。

とりあえず、USE DFLIB となっている行を全部つぶすべし。
そうして出てきた未定義ルーチン・エラーをひとつひとつ調べるべし。




111 名前:107 mailto:sage [2006/12/14(木) 23:36:52 ]
>>110
やはりそうですかー。気が遠くなりそうですが、やってみます。
ありがとうがざいました。

112 名前:デフォルトの名無しさん mailto:sage [2006/12/15(金) 00:48:58 ]
>>111
あぁ、もしLinux上のIntelFortranを使っているなら、DFORTをIFORTに変えて試してみるといい。
今のWin上のIntelfortranではDFORTはIFORTになっている。
Intelが移植用に互換モジュールを用意してくれているかもしれない。

とりあえず一応マニュアルを検索してみるといい。

どっちにしろコマンドライン引数を取るために使っている程度なら、g77でも移植可能だ。


113 名前:102 [2006/12/15(金) 01:26:40 ]
このプログラムだとわざわざDOループのdaの値を何回も変えてよりP3の値が3に近い値のda
を求めなければいけないのですが、それを実行一回で求めたいのですがどうしたら
よいでしょうか?

114 名前:デフォルトの名無しさん mailto:sage [2006/12/15(金) 01:55:18 ]
>>113
見た感じ非線形だから1回は無理で、はさみうち法とかでイターレションするしかナインでないの?
というかもし1回で求まるなら、はじめっからdaについての式を逆関数的に定式化してあるのでは?


115 名前:デフォルトの名無しさん mailto:sage [2006/12/15(金) 10:24:04 ]
>>112
なるほど。試してみます。
ありがとうございました!!

116 名前:デフォルトの名無しさん [2006/12/17(日) 09:35:52 ]
下記のプログラムについて、fortran77で行うと実行できないのです><
何がおかしいですかね><
よければ、教えていただけませんか?


dimension ia(100),vw(12),ivw(12),vw2(100)
data am,ix,n/1.0 ,0 ,100/

call poiss(am,ix,ia,n,vw,ivw,vw2)

write(*,*) (ia(i),i=1,n)

stop
end
 
  subroutine poiss(am,ix,ia,n,vw,ivw,vw2)
dimension ia(*),vw(*),ivw(*),vw2(*)
real*8 x,y
data bm/-1.0/


if (am.ne.bm) then
bm=am
vw(1)=exp(-am)
vw(2)=vw(1)*am
m=2.0*am+10.0

do 10 k=2,m-1
vw(k+1)=vw(k)*am/k
10 continue



117 名前:デフォルトの名無しさん [2006/12/17(日) 09:38:41 ]
do 30 k=1,m-1
if (vw(k).ne.0.0) then
x=vw(k)
do 20 j=k+1,m
y=x+vw(j)
if (y.eq.x) go to 50
x=y
vw(j)=x
20 continue
go to 40
end if
30 continue

k=m
40 j=m

50 do 70 i=j,m
vw(i)=1.0
continue

mm=j-k+1
ai=mm






118 名前:デフォルトの名無しさん [2006/12/17(日) 09:39:33 ]
do 70 i=1,mm
af=(i-1)/ai
do 60 while(af.gt.vw(k))
k=k+1
60 continue
ivw(i)=k
70 continue
end if

call unfm2(ix,vw2,n)

do 90 i=1,n
k=vw2(i)*ai+1.0
if(k.eq.mm+1) k=mm
j=ivw(k)
do 80 while(vw2(i).gt.vw(j))
j=j+1
80 continue
ia(i)=j-1
90 continue

return
end


119 名前:デフォルトの名無しさん [2006/12/17(日) 09:40:24 ]
subroutine unfm2(ix,vw2,n)
implicit real*8 (a-h,o-z)
dimension vw2(*)
real*8 aa/32771.d0/,b/1234567891.d0/,x
real*8 c/2147483648.d0/,ci/4.6566128730773925d-10/

x=ix

do 10 i=1,n
x=dmod(aa*x+b,c)
vw2(i)=x*ci
10 continue
ix=x

return
end



120 名前:デフォルトの名無しさん [2006/12/17(日) 12:11:06 ]
>>116
質問の仕方がなっちょらん。
本文からコンパイルエラーなのか実行時エラーなのかすら読み取れない。
エラーメッセージを出すくらいの知恵も無いのか!www 

と叱られるよ




121 名前:デフォルトの名無しさん mailto:sage [2006/12/17(日) 14:34:44 ]
>>116
一応DEGITAL Visual Fortranでコンパイルしてみた。
do 70のループが原因でコンパイルが通らない。
後ろの70を170にして、その前のcontinueしかなかったところを70にした。
コンパイルはできるが実行中にArray Bounds Exceededが起こる。

プログラムの中身を追いかけるのは面倒だからやらない。

122 名前:デフォルトの名無しさん [2006/12/21(木) 18:06:44 ]
高校の課題なんですが、妹に教えられなくて兄としての威厳がかかっていますw
正直、プログラムはさっぱりなんですが、シスアド初級持っているために安請け合いしてしまったのが間違いですたw

以下の問いなんですがよろしくお願いします。


乱数(0〜1までの範囲)を7個発生させなさい。
その発生させた7個の乱数はそれぞれa,b,c,d,e,f,gと定義する。
このa〜gの中で、最大値と最小値をとるものを選別するプログラムを作りなさい。

というものなんですが、例として何か提示していただけませんか?
それを本みながら自分で理解して、妹に説明できるようにしたいんで。

123 名前:デフォルトの名無しさん [2006/12/21(木) 20:01:24 ]
program main
a=rand()
b=rand()
c=rand()
d=rand()
e=rand()
f=rand()
g=rand()
write(6,*) 'MAX ',max(a,b,c,d,e,f,g)
write(6,*) 'MIN ',min(a,b,c,d,e,f,g)
stop
end


124 名前:122 [2006/12/21(木) 21:05:49 ]
>>123

この

write(6,*) 'MAX ',max(a,b,c,d,e,f,g)
write(6,*) 'MIN ',min(a,b,c,d,e,f,g)

の、

maxとminっていうのは最初に宣言しなくても機能するんですか?

125 名前:デフォルトの名無しさん mailto:sage [2006/12/21(木) 21:34:43 ]
>>124
機能する。 最大値/最小値を返す関数の総称名として規格で規定されている。
RAND() の方は使えない可能性がある。

126 名前:122 [2006/12/21(木) 23:31:03 ]
>>125

ありがとうございました。

127 名前:デフォルトの名無しさん mailto:sage [2006/12/22(金) 00:35:10 ]
>>122 Fortran90風には国府田。

PROGRAM chinpoko
IMPLICIT NONE
REAL :: a(7)
CALL RANDOM_NUMBER(a)
WRITE(*, *) 'MAX=', MAXVAL(a)
WRITE(*, *) 'MIN=', MINVAL(a)
STOP
END PROGRAM chinpoko

128 名前:デフォルトの名無しさん [2006/12/23(土) 00:01:40 ]
ちょっと122番と若干かぶる要素があるんだけど、

0から1までの乱数@〜Dまで存在すると仮定する。

問1.乱数をそれぞれ、@*A、@*B、@*C、@*D…C*Dとすべてのパターンの積を求める問題。

問2.問1で得られた結果に100を掛けて、小数点以下を切り捨てる。

問3.問2で得られた全ての結果に対して50より小さい場合残し、そうでない場合は残さない。残すもの、残さないものを分けろ。

問4.問3で得られた結果(残った数値)を以下のような形でまとめるようにせよ。

yi(y1,y2,y3…yn)

問5.問4で得られた結果からその中で最大値を求めよ。


という5つの問題がある。元は一個の問題なんだが、問題文が無駄に長いので俺が解釈して5つに分割してみた。


質問1.問1を俺なりにやってみたんだが、いまいちプログラムがわからないので自力で全て値を@*Aというように入力していったのだが、もっといい方法はないかな?
質問2.50以下のものとそうでないものを分けるとき、使う分は「if」でいいのか?
質問3.yi(y1,y2,y3…yn)とあるけど、これはdimension?使えばいいのか?


わかりにくいけど質問の回答よろしくお願いします。

129 名前:128 [2006/12/23(土) 00:03:53 ]
できたら、問ごとにプログラムの例文作ってくださると理解しやすいのでありがたいです。

130 名前:デフォルトの名無しさん mailto:sage [2006/12/23(土) 01:12:30 ]
>>129
FORTRAN77かFortran90か指定はあるか?




131 名前:129 [2006/12/23(土) 16:36:07 ]
特にないです。今使ってるのは、visual fortranってやつなんでたぶん新しい方だと思う。

132 名前:デフォルトの名無しさん mailto:sage [2006/12/23(土) 18:27:45 ]
>>128
PROGRAM q1
IMPLICIT NONE
INTEGER, PARAMETER :: n = 5
INTEGER :: i, j, k
REAL :: x(n), y(n * (n - 1) / 2)
CALL RANDOM_NUMBER(x)
k = 0
DO i = 1, n
DO j = i + 1, n
k = k + 1
y(k) = x(i) * x(j)
END DO
END DO
PRINT *, x
PRINT *
PRINT *, y
STOP
END PROGRAM q1

大文字はFortranに存在する命令。小文字は自分の定義する変数名など。
RANDOM_NUMBER()は0以上1未満の数を返す。
>>132 インデント

133 名前:デフォルトの名無しさん mailto:sage [2006/12/23(土) 19:13:04 ]
>>128
PROGRAM q2
IMPLICIT NONE
INTEGER, PARAMETER :: n = 5
INTEGER :: i, j, k
REAL :: x(n), w(n * (n - 1) / 2)
REAL, ALLOCATABLE :: y(:)
CALL RANDOM_NUMBER(x)
k = 0
DO i = 1, n
DO j = i + 1, n
k = k + 1
w(k) = x(i) * x(j)
END DO
END DO
w = AINT(w * 100.0)
PRINT *, x
PRINT *
PRINT *, w
STOP
END PROGRAM q2

問2
INT()は切り捨て整数化。AINTは切り捨て整数化を実数で返す。
>>133


134 名前:デフォルトの名無しさん mailto:sage [2006/12/23(土) 22:44:21 ]
PROGRAM q3
IMPLICIT NONE
INTEGER, PARAMETER :: n = 5
INTEGER :: i, j, k
REAL :: x(n), w(n * (n - 1) / 2), y((n * (n - 1) / 2))
CALL RANDOM_NUMBER(x)
k = 0
DO i = 1, n
DO j = i + 1, n
k = k + 1
w(k) = x(i) * x(j)
END DO
END DO
w = AINT(w * 100.0)
k = 0
DO i = 1, n * (n - 1) / 2
IF (w(i) < 50.0) THEN
k = k + 1
y(k) = w(i)
END IF
END DO
print *, 'no. of less than 50 =', k
PRINT *, y(1:k)
STOP
END PROGRAM q3

135 名前:デフォルトの名無しさん [2006/12/24(日) 19:37:56 ]
すいません。FORTRANについての質問かどうかわからないのですが…。
2次元配列DDDD(512,512)に物理量、1次元配列X(512),Y(512)に座標値が
格納されている場合で、Text column形式で出力するならば
OPEN(10,FILE='TXTCOL.D')
DO L = 1, 512
DO K = 1, 512
WRITE(10,*) X(K),Y(L),DDDD(K,L)
END DO
END DO
というのは理解しているのですが、Binary Columns形式ですとどのように
書けば良いのでしょうか?
Binary Matrix形式ですと
OPEN(20,FILE='BINMAT.D',ACCESS='direct',FORM='unformatted',recl=8*512*512)
WRITE(20,rec=1) ((DDDD(K,L),K=1,512),L=1,512)
とすれば良いと勝手に思っています(※これだと座標値が入ってません)。
そもそもバイナリ書き出しを全然、理解しておりません。
どなたかお教えいただけないでしょうか?

136 名前:デフォルトの名無しさん mailto:sage [2006/12/24(日) 20:01:14 ]
>>135
いまいち何をしようとしているのか理解できないが、単にバイナリー形式で配列を
書き出したいなら書式なし形式で出せばいい。
WRITE(10) DDDD
でおk
WRITE(10) X,Y,DDD
なら、Xのすべて、Yのすべて、DDDのすべて
の順で書き出される。

直接アクセス方式は、要素をランダムにアクセスしたいときの形式。
OPEN文のRECLは単一要素の大きさなので512*512はいらない。
あと、処理系によってバイトかワードか単位要素か単位が違うのでマニュアルで
確かめる必要がある。
WRITE文でのRECは位置を指定するので、要素ごとに1づつ増やしてゆかねばならない。

直接アクセス方式はしばらく使ってないから細部を忘れたw
マニュアルを穴があくまでにらんで読んでくれ。

137 名前:デフォルトの名無しさん [2006/12/28(木) 04:22:27 ]
fortran77を使用して、cos(pi*x)を初期値とした熱方程式の時間発展を計算したいと思います。
範囲は0≦x≦1です。
これで作られたdatファイル('C:\out.dat')が(exp((-(pi)**2)*0.05))*cos(pi*x)
とぴったり重なる様な結果を得たいのですが、ずれてしまいます。
改善点をご指摘いただければ幸いです。よろしくお願いします。
program heat equation
implicit none
integer i,j,nstep,n,ndim
parameter (ndim=5001)
double precision flam,f,h,u,c1,c2,fi,t,tmax,dt
double precision a,b,c
dimension a(ndim),b(ndim),c(ndim),u(ndim)
n=10
dt=0.005d0
tmax=0.05d0
h=1.0d0/float(n)
flam=dt/(h**2.0d0)
c1=1.0d0-2.0d0*flam
c2=flam
do 1 i=1,n-1
a(i)=c2
b(i)=c1
c(i)=c2
1 continue
do 2 i=1,n+1
fi=i-1
u(i)=f(fi*h)
2 continue
OPEN(11,FILE='C:\out.dat')
nstep=0.0

138 名前:デフォルトの名無しさん [2006/12/28(木) 04:23:45 ]
3 continue
do 4 i=1,n-1
4 u(i)=b(i)*u(i+1)+a(i)*u(i)+c(i)*u(i+2)
do 5 i=n,1,-1
5 u(i+1)=u(i)
a(1)=0.0d0
b(1)=c1
c(1)=2.0d0*c2
a(n+1)=2.0d0*c2
b(n+1)=c1
c(n+1)=0.0d0
nstep=nstep+1.0
t=dt*float(nstep)
do j=1,n+1
if(mod(nstep,5).eq.0)then
write(11,*) j*h-1.0d0/n,u(j)
end if
end do
write(11,*)' '
if(t+dt.lt.tmax+0.10d0*dt) go to 3
end

function f(x)
implicit none
double precision x,f,pi
pi=3.14159265358979323846264338327950288
f=dcos(pi*x)
end


139 名前:デフォルトの名無しさん mailto:sage [2006/12/28(木) 05:40:21 ]
>>137
ぴったり重なるといっても、どの程度の重なりを求めているのか?
そもそも積分方程式を差分化して数値積分している時点で誤差が溜まってゆくので
相応の誤差はある。見積もり以上の誤差が出ていると言いたいのか?

中身を全く見ていないが、表面だけを見た感じでは
floatは単精度を返すので、DBLEにしてみそとか、
『fi=i-1』も単精度に変換されているだろうといえる。

*昔のコンパイラは倍精度に単精度の数を入れると、あまった桁数にはゴミが入るので
結局数値誤差は単精度と変わりないレベルまで汚染されてしまう。



140 名前:137 [2006/12/28(木) 06:08:19 ]
>>137
ご指摘の点、早速改善させて頂きました。
ありがとうございました。
ずれは、x=1付近ではほとんどありませんがx=0では誤差では済まない程あります。
グラフにすると形もcos関数ではなくなってしまっているのですが、
改善する手立てが見つけられないでいます。



141 名前:デフォルトの名無しさん mailto:sage [2006/12/28(木) 12:22:48 ]
>>137
なんでFORTRAN77のプログラムを1カラム目から書いてるの?
整形めんどい。

142 名前:デフォルトの名無しさん [2006/12/28(木) 13:17:56 ]
>>141
s/^/ /

FORTRANネタじゃないけど、まいっか

143 名前:141 mailto:sage [2006/12/28(木) 14:16:16 ]
どうも。
質問者がこれを暗黙に要求するのはどうかと思いまして。

viを起動して137-138をコピーペーストして、
:%s/^/ /
を実行した後数字のラベルの前の空白を手動で2カラム削除しました。
スレ違いで済みませんが、
「7カラム目が数字の行の先頭を2カラム削除する」というのもsedかviで
自動で処理出来るのでしょうか?

144 名前:デフォルトの名無しさん mailto:sage [2006/12/28(木) 15:03:05 ]
/^......[0-9]/s/^..//

145 名前:141 mailto:sage [2006/12/28(木) 17:33:24 ]
>>144
ありがとうございます。

146 名前:137 [2007/01/01(月) 04:02:41 ]
前述の計算は以下のように書き直したら期待通りのグラフが描けました。
       program heat equation
       implicit none
       integer i,j,nstep,n,ndim
       parameter (ndim=5001)
       double precision flam,f,h,u,c1,c2,fi,t,tmax,dt
       double precision a,b,c
       dimension a(ndim),b(ndim),c(ndim),u(ndim)
       n=10
       dt=0.005d0
       tmax=0.05d0
       h=1.0d0/dble(n)
       flam=dt/(h**2.0d0)
       c1=1.0d0-2.0d0*flam
       c2=flam
       do 1 i=1,n
       a(i)=c2
       b(i)=c1
       c(i)=c2
      1 continue
       do 2 i=1,n+2
       fi=i-2.0d0
       u(i)=f(fi*h)
      2 continue
       OPEN(11,FILE='C:\out.dat')
       nstep=0.0

147 名前:137 [2007/01/01(月) 04:13:19 ]
      3 continue
       do 4 i=1,n
      4 u(i)=b(i)*u(i+1)+a(i)*u(i)+c(i)*u(i+2)
       do 5 i=n+1,1,-1
      5 u(i+1)=u(i)
       a(1)=0.0d0
       b(1)=c1
       c(1)=2.0d0*c2
       a(n+2)=2.0d0*c2
       b(n+2)=c1
       c(n+2)=0.0d0
       nstep=nstep+1.0
do j=1,n+1
if(mod(nstep,5).eq.0)then
write(11,*) j*h-1.0d0/n,u(j+1)
end if
end do
write(11,*)' '
if(t+dt.lt.tmax+0.10d0*dt) go to 3
end
function f(x)
implicit none
double precision x,f,pi
pi=3.14159265358979323846264338327950288
f=dcos(pi*x)
end
しかしcos(pi*x)としていた初期関数を
if(x.le.0.5) f=x
if(x.gt.0.5) f=1.0d0-x
に変えてみると誤った計算結果が出ます。
改善点等みつかりましたらご指摘頂けると幸いです。

148 名前:JAVAはじめました [2007/01/01(月) 07:10:35 ]
初めてプログラミングをやってみましたが本の通りに入力してもできません(涙)。
こんな事聞いたら馬鹿にされるかもしれませんが、これができなければ前に進めません…
ソースコードをコンパイルできません…何度やっても「javacはコマンドとして認識できません」とエラーが起きます。
Java2SDKの本を買って付属のCDROMよりインストールしたのですが、いったいどこが悪いのかわかりません。
ソースコードにはミスないはずなのにPCが壊れてるのかな?…すいませんわかる方教えてください

149 名前:デフォルトの名無しさん mailto:sage [2007/01/01(月) 13:46:06 ]
>>148
このスレッドはFORTRANという(恐らくキミが生まれるより以前から存在する)言語の
スレッドであって、Javaのスレッドではありません。適切なスレッドに映ってください。

一般的に、コンピュータは同じ事を繰り返すと同じ結果が得られます。
何も変更しないで同じ事を繰り返して、↓こういう結果になるのは極めて当然、当たり前。
> 「javacはコマンドとして認識できません」とエラーが起きます。

150 名前:デフォルトの名無しさん [2007/01/05(金) 21:11:49 ]
登山道がAからFまであるとする。

さいころX,Yを二つ振って、出た目の合計でその登山道を選ぶかどうか決定する。

出た目の合計が7以上であるならば、その登山道を登ることができる、とする。

それぞれの結果を、Aは**、Bは**というように表示し、

Aは8だから行ける、Bは5だから行けないと表示するプログラムを作りなさい。


こんなプログラムなんだが、まずサイコロをどう表現したらいいかわからない点、
次に、でた結果を全てAは**、Bは**というように表示すんのはいいんだけどさ、

行くか行かないかを判別するのって一個一個ifで分岐させてやらないかんの?



151 名前:デフォルトの名無しさん [2007/01/05(金) 21:19:44 ]
追加で書かせてもらうと、

if ( AからFの結果 >7) then
write文で表示する
else
データ破棄、


みたいな感じにしたいんだけど…
プログラム全然できなくて日本語はいりまくりですが…
 

152 名前:age [2007/01/05(金) 21:55:11 ]
おお!こんな奇跡のスレがあるとは…
まじヘルプお願いします。


        え-------------か
      /
    い
  /   \
あ       お-------------き
  \
      / 
    う -------------------く


上記のような図がある。
あ〜くはそれぞれ独立した乱数である。

あ〜くで、次の経路の距離は(例:あ−い間はあ・いの積を小数点第3位で四捨五入したもの
)とする。

終点までの全経路の距離の合計を表示し、(例:あ+い+え+か、あ+う+く etc)
もっとも合計距離が小さいものを選ぶプログラムを作りなさい。


なんか上のほうにも似たようなものたくさんあるけど、応用力がなさすぎて作れない…
助けてください…orz

153 名前:デフォルトの名無しさん mailto:sage [2007/01/06(土) 05:53:37 ]
>>152
Fortranで挑戦、途中で挫折・・・
その後、rubyで作成しました。
あいうえおかきく→abcdefghと置き換えて
tree[]にabc...fghの値を入れると答えを見つけてくれます

#------------------------------
a=0; b=1; c=2; d=3; e=4; f=5; g=6; h=7;
null = nil
$ans = []

# a,b,c,d,e,f,g,h
tree=[1,2,3,4,5,6,7,8]

aa = [b,c]; bb = [e,d]
cc = [d,h]; dd = [g,null]
ee = [f,null]; ff = [f,null]
gg = [g,null]; hh = [h,null]

node = [aa,bb,cc,dd,ee,ff,gg,hh]

def len(tree,node,i,n)
tree[i]*tree[node[i][n]]
end

def path(x)
pair = ["a","b","c","d","e","f","g","h"]
pair[x]
end




154 名前:デフォルトの名無しさん mailto:sage [2007/01/06(土) 05:54:24 ]
>>153

def get_ans(node,tree,i,n,sum,_path_)
_next = node[i][n]
if _next == i then
a = [[sum,_path_]]
$ans = $ans + a
elsif _next == nil
nil
else
x = len(tree,node,i,n)
for j in 0..1
get_ans(node,tree,_next,j,sum+x,_path_+path(_next))
end
end
end

# --- main ---
for i in 0..1
get_ans(node,tree,0,i,0,path(0))
end
p $ans.sort
# --- end main ---

#[[27, "ach"], [38, "abdg"], [42, "abef"], [43, "acdg"]]


155 名前:デフォルトの名無しさん mailto:sage [2007/01/06(土) 20:59:33 ]
>>152
real function distance(path)
common /rval/ rval
real rval(10)
integer path(*), i

distance = 1.0
do i=1,1000
if (path(i) .lt. 1) then
distance = int((distance + 0.5) * 100) / 100.0
return
endif
distance = distance * rval(path(i))
enddo
write(*,*) 'internal error'
stop
end


156 名前:デフォルトの名無しさん mailto:sage [2007/01/06(土) 21:00:26 ]
>>155
program f152
common /rval/ rval
real d(6), minimum
integer i, path(10,6)
data path /
+ 1,2,4,6,0,0,0,0,0,0, ! あいえか
+ 1,3,5,2,4,6,0,0,0,0, ! あうおいえか
+ 1,2,5,7,0,0,0,0,0,0, ! あいおき
+ 1,3,5,7,0,0,0,0,0,0, ! あうおき
+ 1,3,8,0,0,0,0,0,0,0, ! あうく
+ 1,2,5,3,8,0,0,0,0,0 ! あいおうく
+ /
real rval(10)
rval(1) = rand(1) ! seed
do i=1,10
rval(i) = rand(0) * 10
enddo
minimum = 3.40282347E+38
do i=1,6
d(i) = distance(path(1,i))
if(minimum .gt. d(i)) minimum = d(i)
write(*,'(i10,f10.3)') i,d(i)
enddo
write(*,'(a10,f10.3)') 'mininum :',minimum
stop
end

157 名前:デフォルトの名無しさん [2007/01/07(日) 17:59:33 ]
質問させてください.
Windowsのgfortranだとコンパイルも実行もできるソースが,
Linux版のgfortranではコンパイルは通るのに実行時に
  At line 91 of file LHDmake.f90
  fortran runtime error;Bad adress
となって止まってしまいます.

LHDmake.f90はMODULEで,その91行目は「WRITE(21,*) ''」という文(データを改行するため)
なのですが,なにがおかしいのでしょうか?

ちなみにUNIT=21は
「OPEN(UNIT=21,file='Result/plotLHD1_polygon.dat',status='unknown')」
です.

コンパイラのバージョンは
Win:gcc version 4.2.0 20060401 (experimental)
Linux:gcc version 4.0.4 20060507 (prerelease) (Debian 4.0.3-3)
です

158 名前:デフォルトの名無しさん mailto:sage [2007/01/07(日) 20:33:38 ]
>>157
こんなのがあった
ttp://gcc.gnu.org/ml/gcc-bugs/2006-01/msg00595.html


159 名前:デフォルトの名無しさん [2007/01/08(月) 16:12:05 ]
「0.360e-10」は0.360×(10の-10乗)ですよね。

eがなくて、「0.360+250」のように表示されるのは
どういう意味でしょうか?

160 名前:デフォルトの名無しさん mailto:sage [2007/01/08(月) 16:51:03 ]
指数が3桁になって表示を切り詰めたんじゃないのかな?



161 名前:デフォルトの名無しさん mailto:sage [2007/01/08(月) 19:14:02 ]
>>159
マニュアルでFORMAT文のところを嫁

162 名前:デフォルトの名無しさん [2007/01/09(火) 15:52:01 ]
すいません…どなたか150の回答お願いします…







[ 続きを読む ] / [ 携帯版 ]

前100 次100 最新50 [ このスレをブックマーク! 携帯に送る ] 2chのread.cgiへ
[+板 最近立ったスレ&熱いスレ一覧 : +板 最近立ったスレ/記者別一覧]( ´∀`)<120KB

read.cgi ver5.27 [feat.BBS2 +1.6] / e.0.2 (02/09/03) / eucaly.net products.
担当:undef