992コメント455KB

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

■ このスレッドは過去ログ倉庫に格納されています

0001デフォルトの名無しさん

2009/01/24(土) 18:32:01

このスレッドは、他のスレッドでは書き込めない超低レベル、
もしくは質問者自身何が何だか分からない質問を勇気を持って書き込むスレッドです。
FORTRAN使いが優しくコメントを返しますが、
お礼はFORTRANの布教と初心者の救済と次期Fortran2008規格でのCOMEFROM文採用をお願いします。

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

●前スレ
くだすれFORTRAN（超初心者用）その3
http://pc11.2ch.net/test/read.cgi/tech/1196384126/

くだすれFORTRAN（超初心者用）その2
http://pc11.2ch.net/test/read.cgi/tech/1164121236/

くだすれFORTRAN（超初心者用）
http://pc8.2ch.net/test/read.cgi/tech/1138063703/

●関連スレ
FORTRAN Ⅳ
http://pc11.2ch.net/test/read.cgi/tech/1163319215/

0266259

2009/06/09(火) 16:56:11

>>262さん　 >>263さん
どうやら、説明が足りなかったようで、すいません。
今日授業に出たところ、
program
real ::r,a
integer ::io
print *,'CTRL-Dを入力すると終了'
do while (io==0)
if (io==0) then
read (5,*,iostat=io) r
a=4*3.141592*r**2
print *,'半径=',r,'面積=',a
end if
end do
stop
end program
が、答えだそうです。
でもこれを実行すると、CTRL-Dで終了するとき、
前回の計算結果が表示されてしまいます。
またまたすいませんが、何を付け加えればいいのでしょうか？

0267デフォルトの名無しさん

2009/06/09(火) 17:09:31

なんかアホらしい授業だけど。
program
real ::r,a
integer ::io
print *,'CTRL-Dを入力すると終了'
do while (io==0)
if (io==0) then
read (5,*,iostat=io) r
if (io==0) then
a=4*3.141592*r**2
print *,'半径=',r,'面積=',a
end if
end if
end do
stop
end program
どう考えてもexitを使うほうが素直だ。

0268デフォルトの名無しさん

2009/06/09(火) 17:11:10

ありゃ、間違えた。
program
real ::r,a
integer ::io
print *,'CTRL-Dを入力すると終了'
do while (io==0)
read (5,*,iostat=io) r
if (io==0) then
a=4*3.141592*r**2
print *,'半径=',r,'面積=',a
end if
end do
stop
end program

readをifの外側に出さないといかん。

0269デフォルトの名無しさん

2009/06/09(火) 17:20:58

ioの初期化も足りないし。ぐだぐだ。

0270デフォルトの名無しさん

2009/06/09(火) 20:24:17

exitの偉大さを実感すると言う狙いなのだろう。

0271デフォルトの名無しさん

2009/06/09(火) 20:31:01

蛇足だが >>259 の
do while ( .true. )
は
do
だけでおｋ

まあ、
do while (io==0)
に直させるためにワザと do while にしたのかもしれないが

0272デフォルトの名無しさん

2009/06/09(火) 21:52:10

do while (io==0)
if (io==0) then
たまにこういうセンスなさすぎなコード書く人いるよね
他人が書いたコードをメンテするはめになったときに、こういうの見ると何とも言えない溜息がでるわ・・・
たとえバグはないとしても

0273デフォルトの名無しさん

2009/06/09(火) 23:17:11

>>272
まぁ待てｗ
Fortran で WHILE を使おうとするだけでも、きわめて開明的な人物ではないか！

0274デフォルトの名無しさん

2009/06/09(火) 23:27:58

> if (io==0) then
これはどこがセンスなさすぎなん？

0275デフォルトの名無しさん

2009/06/10(水) 03:23:56

>>272 は「センス」って言葉を使いたかっただけなんだよ。

0276デフォルトの名無しさん

2009/06/10(水) 11:45:30

>>274
ループの脱出（継続）分岐条件である　io==0　が
WHILE と　IF の二箇所で判断されていることを言ってるんだろ。

>>272の態度はむかつかせるだけから、くだすれ的には反省すべき。

0277デフォルトの名無しさん

2009/06/10(水) 15:02:05

将来コードをいじったときにそのループをスキップする場合があることを考えて
do while (io .EQ. 0)
みたいなのは意図的に書くことはあるけどね。
というか、>>266 のは変数 io がゼロに初期化されることを
前提にしている時点であまりいい回答では無いように思うけどな。

0278259

2009/06/10(水) 16:33:09

>>268さん、>>271さん
どうも、ありがとうございました。
授業はどうもいまいちわかりません(笑）
これからも、ちょくちょく質問させていただきますので、
くだすれの皆さん、よろしくお願いします。

0279デフォルトの名無しさん

2009/06/11(木) 23:43:42

Visual Fortranについて質問してもいいですか？（VS2008スレよりもこちらのほうが適切っぽいので・・・）

VBやC#でVS2005を使ってたのですが、Intel コンパイラv11を買って
Visual Fortranを入れてみたところ、Intellisenseが効きません。
Visual Fortranでは使えないということですか？
これだとVSをFortran開発環境とする利点がかなり消えてしまうと思うのですが・・・

0280デフォルトの名無しさん

2009/06/12(金) 00:01:40

>>279
使えない。
そういうものがあることを知らなければ気にならないから、いますぐ忘れるんだ！

ただしPolyHedoronがLaheyFortran用のインテリセンスっぽい奴を売っている。
ttp://www.lahey.com/LFVSIntegration.htm

0281デフォルトの名無しさん

2009/06/12(金) 00:46:40

>>280
そうですか・・・
計算はLinuxクラスタで走らすけど、開発は使い慣れたWindowsのVSで
っていう人を見たことがあって、それいいなぁと思って入れてみたんですが・・・

Laheyっていうの初めて知りました。intellisenseっぽいポップアップが出てますね。
これやNAGなど他の開発環境で良いのがないか当たってみようかなぁ

ちなみにNAGの↓の比較図ですが、他社製ってやっぱVisual Fortranですかね・・・（苦笑
ttp://www.nag-j.co.jp/nagfor/img/f95comp.jpg

0282デフォルトの名無しさん

2009/06/12(金) 00:57:25

>>281
PGI Visual Fortran というのもあって、これもVSと統合されているのだが、
確か今は米国内の研究機関に所属していればただで使える。
アメリカの研究機関に知り合いがいれば（略

SalfordというかSilverFrostのFTN95もVSと統合できる。
しかし、PGIもSalfordもインテリセンスは無かったような気がする。
最新状況はしらねｗ

IntelFortranは今月中にF2003対応度が高まったv11.1で出るらしいので、
あわてずゆっくりしていってね。

0283デフォルトの名無しさん

2009/06/12(金) 03:33:59

>>282
色々ありがとうございます。

↓を見みると、SalfordのVS統合はintellisenseみたいなのがあるようです。
http://www.silverfrost.com/17/ftn95/ftn95_visual_studio_dotnet_integration.aspx
というかこれ無料っぽくないですか？太っ腹すぎる・・・

PGIも、ざっくり見た感じIntel VFよりは機能が多そうですね。
試用できるようなので試してみます。

しかしIntel VFはせめてキーワード補完くらいあればいいのに・・・
性能第一で見やすさ・作りやすさは二の次ってことなんでしょうか（ある意味硬派？

0284279

2009/06/12(金) 21:26:08

さっそく以下の3つを試しました。

Silverfrost FTN95 Express
PGI Visual Fortran
NAG Fortran Builder

SilverfrostとPGIはVS組み込みタイプで、IVFと同じくintellisense等の高度なものはありませんでした。
それ以外も、パッと見特に違いはなさそうで、結局は好みの問題ということになるんでしょうか。
NAGのIDEはVS組み込みではなくオリジナルで、唯一これだけキーワード補完ができましたが、
変数や関数の補完は×でした。

うーん、こうなるとこのままIVFを使い続けたほうがよさそうです・・・

0285279

2009/06/12(金) 21:44:13

ググっていたところ、IntelのForumでどんぴしゃなスレを見つけました。
http://software.intel.com/en-us/forums/intel-visual-fortran-compiler-for-windows/topic/57076/

皆思うことは同じのようですねｗ
Steveおじさんも要望をよく理解してくれてるようですが、これまではVS2008への対応で
手一杯だったようです。MSが仕様をコロコロ変えることも妨げの原因だとぼやいてました。

結局、何が根本の原因なのかははっきりと書かれてませんでしたが、IVF以外のVS組込みIDEも
IVFと似たような現状なのを考えると、VS側のほうでAPIみたいな道具を提供してくれないと
Intellisenseのような高度機能は実装できない、ということなのかもしれません。

Steveさんが言うに、C++,C#,VBのVSで備わっていたこれらの機能が「Visual Studio Shell」
で分離されたので、Fortran上にも実装できそう、とのことです。
早くそうなって欲しいところですが、先月のスレへのSteveさんのレスを見る限り、もう少しかかりそうですね・・・
http://software.intel.com/en-us/forums/intel-visual-fortran-compiler-for-windows/topic/65546/

0286デフォルトの名無しさん

2009/06/13(土) 16:09:35

確かにfortranerの開発って時代遅れな感があるな

0287デフォルトの名無しさん

2009/06/14(日) 22:04:16

>>285
スティーブって森の熊さんホモ動画の熊さんの一匹に似てるよね。

0288デフォルトの名無しさん

2009/06/15(月) 22:36:44

Fortranの開発環境のヘボさにうんざりしてCに逃げる人も多いんだろうな
ただ、今後のトレンドを考えるとやはりCは…

0289デフォルトの名無しさん

2009/06/15(月) 22:46:01

紙テープ・パンチカード・ラインエディター上がりも多いＦｏｒｔｒａｎ界では、スクリーンエディターで
シンタックス・カラーリングがあるだけでも極楽浄土のけだるい怠けた世界に見えるのさ。

0290デフォルトの名無しさん

2009/06/15(月) 23:20:51

これから並列数は増える一方だから、最適化は少しの妥協も許されなくなるんだよな
そうなるとFortran一択しかなくなる…
MSがVisualFortran出してくれれば(^q^)

0291デフォルトの名無しさん

2009/06/15(月) 23:30:16

>>290
っ [MS Fortran PowerStation]

0292デフォルトの名無しさん

2009/06/16(火) 00:33:41

>>291
そう、それをもう一度頼む！(><)
MSのFortran部隊はどこ行っちまったんだ？

0293デフォルトの名無しさん

2009/06/16(火) 00:49:49

>>292
Netscapeがブラウザー戦争を仕掛けたせいで、VC++部隊以外が全部IE開発に
回されたという噂が当時あった。

まぁMSの歴史をたどると、BASICの次に作って売ったのはFORTRANだから
MS Visual FORTRANは原点回帰といえなくもないが、VBすら捨ててしまうMS様は
容赦ないのだｗ

でも、高い金出してコンパイラ律儀に買ってくれるのはFORTRANユーザーくらい？

0294デフォルトの名無しさん

2009/06/16(火) 00:55:45

>>292
MS → DEC → Compaq → Intel

0295デフォルトの名無しさん

2009/06/16(火) 13:53:00

>>294
MSからDECに流れたとはしらなんだ
IntelとMSがもっとみっちり連携してくれればいいんだけど、またWintelって非難が起こりそうだな

0296arrayの地獄

2009/06/16(火) 22:24:17

時間がかかって仕方ないのだが何とかなりませんか？
ifort 9.1です。

real(8)::x(25,25),y(50,50)
x=なんたら
y=0
y(1:25,1:25)=x(1:25,1:25)
y=cshift(y,12,1)
y=cshift(y,12,2)

特にy=xのところ。
なんでこんなんに時間食う？？
（メモリ書き込みに難儀してるのはわかるが、
y全体がキャッシュ上に乗るから、キャッシュで整えてから書けば時間かからんはず・・）

0297デフォルトの名無しさん

2009/06/16(火) 22:27:38

>>296
それだけじゃ分からんっつの
時間食うってどれくらい？

0298arrayの地獄

2009/06/16(火) 22:46:24

９２ｘ４９回で数秒のオーダー
＠xeon5160？ 1core使用
なんかおかしいな・・もうちょっと調べます。

0299デフォルトの名無しさん

2009/06/16(火) 22:49:48

ド初心者ですが質問です。コンパイルできたファイルを実行すると
fmt: end of file
apparent state: unit 10 named inp-2.dat
last format: (F9.5)
lately reading sequential formatted external IO
となるのですが何が問題なんでしょうか？プログラムは外部データ（inp-2.dat）から数値を読み最大値を求めるものです
数値は実数です

0300デフォルトの名無しさん

2009/06/16(火) 22:56:05

>>296
>y(1:25,1:25)=x(1:25,1:25)

なぜ素直にｙ＝ｘとしない。　

多分
ALLOCATE（TMP(26,26))
TMP=x
ｙ＝TMP
DEALLOCATE（TMP)
のように展開されるから遅いんだろ。

0301デフォルトの名無しさん

2009/06/16(火) 22:59:20

>>300
TMP(25,25)だったｗ　0:25と勘違いしたよ。1:25とか俺的にはあり得ないｗ

>>299
データ数が足りないように見受けられる。
F9.5に当たるところでデータが足りないか、変数の数が合わないかしてるのでは？

0302296です。

2009/06/16(火) 23:06:23

[hoge@xeon prog20]$ cat test.f90
program test
implicit none
integer,parameter::nx=25,ny=50
real(8)::x(nx,nx),y(ny,ny)
integer::i
x=reshape((/1:25*25/),shape(x))
do i=1,100000
call sub1(x,y,nx,ny)
call sub2(x,y,nx,ny)
end do
end program test
subroutine sub1(x,y,nx,ny)
implicit none
real(8)::x(nx,nx),y(ny,ny)
integer::nx,ny
y=0
y(1:nx,1:nx)=x(1:nx,1:nx)
y=cshift(y,12,1)
y=cshift(y,12,2)
end subroutine
subroutine sub2(x,y,nx,ny)
implicit none
real(8)::x(nx,nx),y(nx,nx) !nyでなくnx
integer::nx,ny
y=0
y(1:nx,1:nx)=x(1:nx,1:nx)
y=cshift(y,12,1)
y=cshift(y,12,2)
end subroutine　続く

0303デフォルトの名無しさん

2009/06/16(火) 23:07:10

>>301
データ数ですか。それは入力元のファイルの値に問題ありということですか？プログラムは
IMPLICIT REAL(A-H,O-Z)
DIMENSION X(100)
OPEN(UNIT=10,FILE='inp-2.dat',STATUS='old',FORM='FORMATTED')
OPEN(UNIT=11,FILE='list-2.dat',STATUS='UNKNOWN',FORM='FORMATTED')
READ(10,100)N
DO 20 I=1,N
READ(10,200)X(I)
20 CONTINUE
IMAX=0.0
DO 30 I=1,N
IF(X(I).GT.IMAX)MAXN=I
IF(X(I).GT.IMAX)IMAX=X(I)
30 CONTINUE
WRITE(11,300)'MAX=',IMAX,'NO.=',MAXN
100 FORMAT(A5)
200 FORMAT(F9.5)
300 FORMAT(A4,F10.5,A4,I3)
CLOSE(10)
CLOSE(11)
STOP
END
となってます。サンプルプログラムを改変しただけなのでおかしいところが多いと思います。入力ファイルはとりあえず実数を5つ縦に並べてあります

0304296です。

2009/06/16(火) 23:07:11

[hoge@xeon prog20]$ ifort test.f90 -p
[hoge@xeon prog20]$ time a.out
real 0m2.359s
user 0m2.357s
sys 0m0.002s
[hoge@xeon prog20]$ gprof -b
Each sample counts as 0.01 seconds.
% cumulative self self total
time seconds seconds calls s/call s/call name
71.83 1.53 1.53 100000 0.00 0.00 sub1_
28.17 2.13 0.60 100000 0.00 0.00 sub2_
0.00 2.13 0.00 1 0.00 2.13 MAIN__

0305296です。

2009/06/16(火) 23:13:31

みにくてすみませんでした。

CPUはcore2の6600 2.4Gでcashe＝4Mでした。

yとxの寸法が違うので、y=xにできないはずです。

逆フーリエ変換にかけるとき外を0にしておきたいので、
こんなになっています。
実プログラムでは逆フーリエ変換にかかる時間と
同じくらいの時間がかかって大損しています。

0306296です。

2009/06/16(火) 23:18:31

まとめると、
y(1:25,1:25)=x(1:25,1:25)
を100000回実行すると、
real(8)::x(25,25),y(50,50) のとき　=>1.53秒かかる。
real(8)::x(25,25),y(25,25) のとき　=>0.60秒かかる。
2.5倍もかからんでええやろ
ということです。

0307296です。

2009/06/16(火) 23:27:05

>>299
私もできる範囲で手伝い。
100 FORMAT(A5) が気になります。
実際に実行させたらNにどんな値が入るでしょうか？
思っている値かどうかwrite文で確かめるのが吉かと思います。

0308デフォルトの名無しさん

2009/06/16(火) 23:32:44

>>307
このままだと実行すらできないんですよね・・・　
READ（～）N　　のあとに
N=5　　としてみましたが実行エラーの文は変わりませんでした

0309デフォルトの名無しさん

2009/06/16(火) 23:40:27

>>296
corei7のifort 11ですが

$ cat test.f90
integer,parameter::n=25
real(8)::x(n,n),y(2*n,2*n)
x=0
y=1
!y=x
y(1:n,1:n)=x(1:n,1:n)
end

$ time ./a.out
real 0m0.001s
user 0m0.000s
sys 0m0.001s

n=25000だと
real 0m1.974s
user 0m0.754s
sys 0m1.217s

何か環境おかしいんじゃないんでしょうか？他のホストで試したどうっすかね

0310296です。

2009/06/16(火) 23:44:15

>>299
適当にinp-2.datに
5
1
12
123
1234
12345
123456
と入れてます。
read(～)Nのあとにwrite(*,*)Nとすれば、ここまでは走りませんか？

これで走らせると、私の場合はN= 538976288と、とんでもないことになってました。
数字を数として読み込みたいのに、文字として読み込んだらまずいですよね、やっぱり。
A5はI5ぐらいかと思いますが、どうでしょう？

0311296です。

2009/06/16(火) 23:49:47

>>309
情報ありがとうです。
corei7＋ifort 11ですか・・
なるほど。別のマシンで試してみます。

0312デフォルトの名無しさん

2009/06/16(火) 23:54:28

>>310
inp-2.datに小数点以下１桁を含んだ実数を適当に5つ入れて、writeしてみたところN＝1073741824となってました。
フォーマットをi5にすると実数を受け取れないのでF5.5にしてみたのですが、Nは上記の通りです。何が起こってるんですかね・・・？

0313296です。

2009/06/16(火) 23:55:32

>>309
すんません。メインルーチンで
do i=1,100000
と回してますので、timeの結果が直接比較できないです。
すみませんが追加でテストお願いできませんか？

0314デフォルトの名無しさん

2009/06/17(水) 00:14:43

>>313
いやメインルーチンとかじゃなくて
>>309でテストすりゃいいじゃん

0315296です。

2009/06/17(水) 00:18:42

>>299
Nはデータの個数ですよね。実数じゃまずいのでは？
inp-2.datは
データの個数
データ１
データ２
：
というような構造じゃありませんか？
少なくとも最初のデータだけは整数かな。
それとIMAXを実数宣言しなきゃいけないようです。

0316デフォルトの名無しさん

2009/06/17(水) 00:27:14

>>315
datファイルにデータの個数を入れてなかったんですね・・・やっとわかりました！ありがとうございます。実行できました
プログラム自体はうまくいってないんですが・・・ｗ
実行結果がMAX＝０．０００・・・、No.＝２となっているんでうまくいってないですね

0317296です。

2009/06/17(水) 00:33:29

>>314
お世話になってます。
確かにこちらでも結果0秒になりました。

いやしかし、たぶん0.000,0060秒（6マイクロ秒）くらいかかっているのが、
timeコマンドの分解能不足で測れてないだけの様な気がしてますです。
実際は100000ｘ100000回とかもっとぐらい実行するので、何とかしたいです。

>>316
IMAXの方もよろしく。

0318309

2009/06/17(水) 00:39:12

>>317
問題点を整理（できればテストコードを用意）して頂ければ協力できますよ
他のマシンもいくつかあります

0319デフォルトの名無しさん

2009/06/17(水) 00:53:40

>>317
実数宣言とは何でしょうか？imax=0.0ではダメなんですよね？

0320300

2009/06/17(水) 01:02:12

>>305
>yとxの寸法が違うので、y=xにできないはずです。
すまんこｗ　ようつべ見ながらだったので、よく見なかったｗｗ

こんなものに時間がかかるも糞もないだろと思ったが、
１０＾５＊１０＾５回も繰り返すなら気になるかｗ

部分配列でメモリーが連続していないので、コンパイラのバージョンが古いと
余計なテンポラリを暗黙にとっている可能性が高いと思う。

DO i = 1, 25
y(1:25, i) = x(:, i)
END DO
と、でもすれば早くなるべ。最適化でアンローリングとかもあろうし。

0321296です。

2009/06/17(水) 01:07:49

>>318
ありがとうございます！
テストコードを302に書き込みました。これでお願いできたらありがたいです！

問題点のほうは、まとめると、
302のプログラムで
y(1:25,1:25)=x(1:25,1:25)
を100000回実行すると、
real(8)::x(25,25),y(50,50) のとき　=>1.53秒かかる。
real(8)::x(25,25),y(25,25) のとき　=>0.60秒かかる。
「2.5倍もかからんでええやろ」
ということです。

おわかりいただけたでしょうか？

ifort10.1を試みてみました。XEON5160です。
% time a10.out　　　　　：ifort10.1
1.288u 0.000s 0:01.31 97.7%
% time a09.out　　　　　：ifort9.1
1.916u 0.004s 0:01.96 97.4%
だいぶ差があるなあ・・というところです。

0322296です。

2009/06/17(水) 01:16:14

>>319

real imax
を宣言文のところに追加です。
FORTRANの慣習ではI～Nから始まる変数名は整数変数に使うことになっているので、
imaxのままだと動くけどちょっと減点かな？
imax->xmaxにするほうがいいかも。この場合real xmaxは不要。

>>320
ご指摘ありがとうですです。
さっそく試してみますです！

0323デフォルトの名無しさん

2009/06/17(水) 01:21:31

>>322
うまくいきました！　なんとなくfortranがわかりはじめてきました！
ありがとうございます

0324デフォルトの名無しさん

2009/06/17(水) 01:30:04

初心者は呪文のようにimplicit noneを宣言しておきましょう。

0325309

2009/06/17(水) 01:39:38

>>321
こんな感じでした。時間はreal timeです

Corei7 920、ifort 11.0
real(8)::x(25,25),y(50,50)・・・0m0.799s
real(8)::x(25,25),y(25,25)・・・0m0.255s

Xeon5472、ifort 10.1
real(8)::x(25,25),y(50,50)・・・0m0.922s
real(8)::x(25,25),y(25,25)・・・0m0.280s

Opteron2218、pgf90 6.2
real(8)::x(25,25),y(50,50)・・・0m9.821s
real(8)::x(25,25),y(25,25)・・・0m4.372s

>>300さんの指摘通りってことですかね
Opteronが糞遅いのが気になりますが・・・（；^ω^）

0326296です。

2009/06/17(水) 01:56:29

>>309様（325）
ありがとうございます！！お手数掛けました。

COREi7万歳というところですか。
320様のご指摘も含めていろいろ検討したいのですが、
今晩はこの辺で失礼します。
（いろいろありがとうございました。また明日頑張ります。）

0327デフォルトの名無しさん

2009/06/17(水) 03:05:28

corei7はかなり評判いいよね
opteronを完全に引き離した

0328デフォルトの名無しさん

2009/06/17(水) 03:11:03

>>325
OpteronというよりPGIのコンパイラのせいのような気がする。
オプションがどうなってるのかわからんし。
最新版は8.0だし。

でもたしかにi7はなかなかいい。

0329309

2009/06/17(水) 13:56:34

>>328
＞OpteronというよりPGIのコンパイラのせいのような気がする。
その通りでしたｗ
>>325のコンパイルは全てオプション指定なしでして、調べたところ
ifortはデフォルトで-O2が、pgf90は-O1がかかるようです。

-O2を付けた場合↓
real(8)::x(25,25),y(50,50)・・・0m1.306s
real(8)::x(25,25),y(25,25)・・・0m0.784s

0330デフォルトの名無しさん

2009/06/17(水) 20:36:22

>>301
＞1:25とか俺的にはあり得ないｗ
Cみたいに0から開始ってこと？
Fortranだと1から始めるのが一般的だと思うけど

0331デフォルトの名無しさん

2009/06/17(水) 23:16:26

>>330
いや、：だけで済ますとか、それすら書かないと言うこと。

0332デフォルトの名無しさん

2009/06/18(木) 18:02:57

>>331
そういうことね

で、実際どれがいいんだろうか。例えば↓の３つでは
① x(1:n) = y(1:n)
② x(:) = y(:)
③ x = y
自分は②が良いと思ってる。①は冗長だし、③は簡潔だけど配列であることが自明じゃないので。
問題はどれが一番速いかなんだけど・・・
↓のサイトでは「並列最適化コンパイラの動作を助け、バグを少なくするため」①が推奨されてるし
ttp://www.mri-jma.go.jp/Project/mrinpd/coderule.html
でも、これって>>300氏の言ってることとは反対だよね。

0333デフォルトの名無しさん

2009/06/18(木) 19:19:25

例えば、配列 data1 から data4 があって、次々にsubroutine に引き渡したい場合ですが、

call hoge(data1)
call hoge(data2)
call hoge(data3)
call hoge(data4)

しないで、

do i=1,4
call hoge( )
end do

の様にまとめる事って加能ですか？　もし可能ならやり方教えてください。

0334デフォルトの名無しさん

2009/06/18(木) 19:22:10

data1 から data4 は既に７次元配列なので、それらをまとめる事はできないです。

0335デフォルトの名無しさん

2009/06/18(木) 19:22:23

>>333
普通は配列をdata(n,4)と2次元にして
do i=1,4
call hoge(data(:,i))
enddo

0336デフォルトの名無しさん

2009/06/18(木) 19:25:49

>>334
7次元配列ってすごいなｗ
Fortranの配列って7が限界だったっけ？それなら無理じゃないかな
あきらめれ

0337デフォルトの名無しさん

2009/06/18(木) 23:01:17

>>332
新しいコンパイラでは、どれを書いても最適化で同じになるらしいが、
初期のFortran90コンパイラでは、n：nnみたいな指定子が入ると、
テンポラリをとることが多かったみたい。

それは多分一般的には、n：nn：mみたいな不連続の番地をとる可能性があるためだと思う。

>>333
ポインターの配列を、TYPEによって定義しておいて、それに配列を割り付けておいて、
それを引数に置けばいいんじゃないかな？

まぁポインターに割り当てるときに、だらだらとdata1~4を並べないと駄目なので
ゴミの位置が移動するだけだがｗ

0338337

2009/06/19(金) 01:27:21

MODULE m_billy
IMPLICIT NONE
INTEGER, PARAMETER :: nmax = 2
TYPE :: t_billy
INTEGER, POINTER :: mm(:, :, :, :, :, :, :)
END TYPE t_billy
CONTAINS
SUBROUTINE anal(mm)
INTEGER, INTENT(OUT) :: mm(:, :, :, :, :, :, :)
INTEGER :: i1, i2, i3, i4, i5, i6, i7
DO i1 = 1, nmax
DO i2 = 1, nmax
DO i3 = 1, nmax
DO i4 = 1, nmax
DO i5 = 1, nmax
DO i6 = 1, nmax
DO i7 = 1, nmax
mm(i1, i2, i3, i4, i5, i6, i7) = i1 + 10 * (i2 + 10 * (i3 + 10 * (i4 + 10 * (i5 + 10 * (i6 + 10 * i7)))))
END DO
END DO
END DO
END DO
END DO
END DO
END DO
RETURN
END SUBROUTINE anal
END MODULE m_billy

0339337

2009/06/19(金) 01:28:03

!====================================================
PROGRAM herrington
USE m_billy
IMPLICIT NONE
INTEGER, TARGET :: data1(nmax, nmax, nmax, nmax, nmax, nmax, nmax), data2(nmax, nmax, nmax, nmax, nmax, nmax, nmax)
INTEGER, TARGET :: data3(nmax, nmax, nmax, nmax, nmax, nmax, nmax), data4(nmax, nmax, nmax, nmax, nmax, nmax, nmax)
TYPE (t_billy) :: kazuya(4)
INTEGER :: i
kazuya(1)%mm => data1
kazuya(2)%mm => data2
kazuya(3)%mm => data3
kazuya(4)%mm => data4
DO i = 1, 4
CALL anal( kazuya(i)%mm )
END DO
print *, data1
STOP
END PROGRAM herrington

0340333

2009/06/19(金) 01:58:08

サンクス。でも結構大掛かりだな。。。

real*8 data1(n,n,n,n,n,n,n), data2(n,n,n,n,n,n,n),.......
char*1 cnum
char*5 target

do i = 1, 4
write( cnum, '(i1)' ) i
target = "data" // cnum
call hoge( target )
end do

みたいにできんもんかね。。。 target と data$ は型が違うがな。

0341337

2009/06/19(金) 02:24:46

>>340
77の掟破りでよければ、

real*8 data1(n,n,n,n,n,n,n), data2(n,n,n,n,n,n,n), data3(n,n,n,n,n,n,n), data4(n,n,n,n,n,n,n)
real*8 tmp(n**7, 4)
equivalence (tmp(1, 1), data1), (tmp(1, 2), data2), (tmp(1, 3), data3), (tmp(1, 4), data4)
do i = 1, 4
call hoge( tmp(i) )
end do

こんな感じでEUIVALENCEでいけるかな？

0342デフォルトの名無しさん

2009/06/19(金) 03:00:32

>>341
call hoge( tmp(1, i) )
の過ち。

ハルヒ新OP糞杉ｗｗｗ

0343デフォルトの名無しさん

2009/06/19(金) 10:22:30

>>334
本当に7次元いるの？　変数設計を再検討した方がいいんじゃない？
座標3次元＋時間で、4次元はざらにあるが。

0344デフォルトの名無しさん

2009/06/19(金) 14:06:33

>>343の言うとおり。

変数設計は、アルゴリズム的な自然さ（人間にとっての見やすさ）と速度の
両方の観点から決めるべき。

0345デフォルトの名無しさん

2009/06/20(土) 20:01:44

7次元配列の数値計算ってどういう分野で出てくるの？

0346デフォルトの名無しさん

2009/06/20(土) 21:51:11

流れぶった切ってすいません。
n×n 次元行列　M と　n 次元ベクトルa の積 b = Ma を計算するプログラムを
作っていたのですが、詰まってしまいました。
どこが、おかしいか教えていただけないでしょうか？
program
integer ::i,k,n
integer,parameter ::limit=5
real ::M(1:limit,1:limit),a(1:limit),b(1:limit)
do
print *,'nの値は？'
read *,n
if ((n<=limit).and.(n>=1)) then
exit
end if
end do
!
print *,'行列Mの成分'
do i=1,n
print *,i,'行目?'
read *,(M(i,k),k=1,n)
end do
print *,'ベクトルaの成分'
do k=1,n
read *,a(k)
end do

0347346

2009/06/20(土) 21:52:41

!>>346の続き
do i=1,n
do k=1,n
b(0)=0
b(k)=b(k)+M(i,k)*a(k)
end do
end do
print *,'b:'
do k=1,n
print *,b(k)
end do
stop
end program

よろしくお願いします

0348346

2009/06/20(土) 21:57:18

たびたびすいません。
上のプログラムはfortran90を使っています。
コンパイルする際は、ifort (プログラム名）.f90 -o (プログラム名）
としています。

0349デフォルトの名無しさん

2009/06/20(土) 21:59:30

>>346
全く見ていないが、b(0)=0　がおかしい。
b（k）=0.0 でかつ一個DOLOOPの外。

0350デフォルトの名無しさん

2009/06/21(日) 10:45:18

ループのところ
do i=1, n
b(i) = 0
do k=1,n
b(i) = b(i) +M(i, k) * a(k)
end do
end do

直接関係ないけど、Mの行と列を入れかえたほうがいいかも。

0351346

2009/06/21(日) 18:17:49

>>350さん
ありがとうございました。
うまくいきました！

0352デフォルトの名無しさん

2009/06/23(火) 11:08:29

PROGRAM GAUSS_SOLVER
IMPLICIT NONE

REAL,DIMENSION(100,100)::A !行列A（最大100元まで解ける）
REAL,DIMENSION(100)::B,X !右辺Bと解Xのベクトル
INTEGER::N !元数
INTEGER::I,J,K !制御変数

WRITE(*,*)'このプログラムは連立一次方程式AX=Bのい解Xを求めます'
WRITE(*,*)'AandB.TXTの入力フォーマットは次のようです'
WRITE(*,*)'たとえば,３元の連立方程式の場合'
WRITE(*,*)'3 (元数を与えてください)'
WRITE(*,*)'A11 A12 A13 B1 (データは空白やタブで別けてください)'
WRITE(*,*)'A21 A22 A23 B2 (AIJとBIには実数を与えてください)'
WRITE(*,*)'A31 A32 A33 B3'
WRITE(*,*)'解XはX.TXTに出力されます'

OPEN(1,FILE='AandB.TXT')
READ(1,*) N
DO I=1,N
READ(1,*) (A(I,J),J=1,N),B(I)
END DO
CLOSE(1)

DO I=1,N-1
DO J=I+1,N
DO K=I+1,N
A(J,K)=A(J,K)-A(J,I)*A(I,K)/A(I,I)
END DO
B(J)=B(J)-B(I)*A(J,I)/A(I,I)
END DO
END DO

0353デフォルトの名無しさん

2009/06/23(火) 11:10:04

DO I=N,1,-1
X(I)=B(I)/A(I,I)
DO J=I+1,N
X(I)=X(I)-A(I,J)*X(J)/A(I,I)
END DO
END DO

OPEN(1,FILE='X.TXT')
DO I=1,N
WRITE(1,*) 'X(',I,')=',X(I)

END DO
CLOSE(1)

STOP
END PROGRAM GAUSS_SOLVER

0354デフォルトの名無しさん

2009/06/23(火) 11:20:09

このプログラムをもとにN次行列の逆行列をもとめるプログラムを作るのですが
どうしたらいいかわかりません。教えてくださると助かります。。

0355デフォルトの名無しさん

2009/06/24(水) 00:59:41

ガウス・ジョルダン消去法でぐぐれ

ちなみに逆行列を求める時は元の配列に逆行列を
入れる事が出来て効率が良い

0356デフォルトの名無しさん

2009/06/24(水) 22:19:20

subroutine gauss_jordan(a,b,flag)
implicit none

integer,parameter::N=3
double precision ::a(N,N),b(N)
integer ipv,i,j
double precision :: inv_pivot,temp
double precision :: big
integer pivot_row,row(N)
integer flag

do ipv=1,N
!ｺﾇﾂ酖ﾍﾃｵｺ・ big=0.0
do i=ipv,N
if(abs(a(ipv,i)) > big) then
big = abs(a(ipv,i))
pivot_row = i
endif
enddo

if(abs(big) < 0.0001) then
write(6,*) "A is singular matrix"
flag = 0
return
endif

0357デフォルトの名無しさん

2009/06/24(水) 22:20:34

row(ipv) = pivot_row

!ｹﾔ､ﾎﾆ､・ﾘ､ｨ
if(ipv .ne. pivot_row) then
do i=1,N
temp = a(i,ipv)
a(i,ipv) = a(i,pivot_row)
a(i,pivot_row) = temp
enddo
temp = b(ipv)
b(ipv) = b(pivot_row)
b(pivot_row) = temp
endif

! ﾂﾐｳﾑﾀｮﾊｬ=1(･ﾔ･ﾜ･ﾃ･ﾈｹﾔ､ﾎｽ靉)
inv_pivot = 1.0/a(ipv,ipv)
a(ipv,ipv) = 1.0

0358デフォルトの名無しさん

2009/06/24(水) 22:21:30

do j=1,N
a(j,ipv) = a(j,ipv) * inv_pivot
enddo

b(ipv) = b(ipv) * inv_pivot

! ･ﾔ･ﾜ･ﾃ･ﾈﾎ・0(･ﾔ･ﾜ･ﾃ･ﾈｹﾔｰﾊｳｰ､ﾎｽ靉)
do i=1,N
if(i .ne. ipv) then
temp = a(ipv,i)
a(ipv,i) = 0.0
do j=1,N
a(j,i) = a(j,i) - temp*a(j,ipv)
enddo
b(i) = b(i) - temp*b(ipv)
endif
enddo

0359デフォルトの名無しさん

2009/06/24(水) 22:23:23

enddo

! ﾎﾎﾆ､・ﾘ､ｨ(ｵﾕｹﾔﾎ・
do j=N,1,-1
if(j .ne. row(j)) then
do i=1,N
temp = a(j,i)
a(j,i) = a(row(j),i)
a(row(j),i) = temp
enddo
endif
enddo

flag = 1
return
end subroutine

!### TEST MAIN ######
! program MAIN
! implicit none
! integer,parameter::N=3
! double precision :: a(N,N),b(N)
! integer :: flag

0360デフォルトの名無しさん

2009/06/24(水) 22:25:41

! a(1,1) = 0.0
! a(2,1) = 1.0
! a(3,1) = 0.0
! a(1,2) = 1.0
! a(2,2) = 0.0
! a(3,2) = 0.0
! a(1,3) = 0.0
! a(2,3) = 0.0
! a(3,3) = 1.0
!
! b(1) = 3.0
! b(2) = 4.0
! b(3) = 5.0

! call gauss_jordan(a,b,flag)
! write(6,*) b(1),b(2),b(3)

! stop
! end program

0361デフォルトの名無しさん

2009/06/24(水) 22:27:01

UTF-8か何かで貼ったのか？
文字化けしてて日本語部分読めないぞ

0362デフォルトの名無しさん

2009/06/24(水) 22:27:21

これのことですか？長々とすみません

0363デフォルトの名無しさん

2009/06/24(水) 22:29:24

ダウンロードしたものをノートパッドで開いてコピペしたんですが。。

0364デフォルトの名無しさん

2009/06/24(水) 22:33:53

で、わかったっしょ？ガウスジョルダン消去法を使うと逆行列が
求められるという事が

もちろん行列式が零だと潰れてしまうが

0365デフォルトの名無しさん

2009/06/24(水) 22:40:14

半分くらい・・・（汗

0366デフォルトの名無しさん

2009/06/24(水) 22:42:15

Cのソースで良いなら

http://nyan11.ciao.jp/B/NA/EWSNA/linkANNA.html

のINV1.Cが逆行列を実際にガウスジョルダン消去法で解くプログラムだ
これを参考にFORTRANに移植するとよい

■ このスレッドは過去ログ倉庫に格納されています