zyj8119

木蟲 (著名寫手)

應(yīng)助: 65 (初中生)
貴賓: 0.003
金幣: 915.1
散金: 1440
紅花: 35
帖子: 2936
在線: 1329.4小時(shí)
蟲號(hào): 664177
注冊(cè): 2008-11-29
性別: GG
專業(yè): 理論和計(jì)算化學(xué)

[交流] 【求助】問一個(gè)FORTRAN程序已有3人參與

CODE:

 subroutine dsrrt(a,xr,xi,n,m,l,b)

        dimension a(n),xr(m),xi(m),b(n)

        if(abs(a(1))+1.0.eq.1.0)then

       l=0

         write(*,5)

         return

        end if

5     format(1x,'err')

      l=1

      k=m

      is=0

        w=1.0

        do 10 i=1,n

10    b(i)=a(i)/a(1)

20    pp=abs(b(k+1))

      if(pp.lt.1.0e-12)then

        xr(k)=0.0

        xi(k)=0.0

        k=k-1

        if(k.eq.1)then

         xr(k)=-b(2)*w/b(1)

         xi(k)=0.0

         return

        end if

        goto 20

        end if

        q=pp**(1.0/k)

        p=q

        w=w*p

        do 30 i=1,k

        b(i+1)=b(i+1)/q

        q=q*p

30    continue

      x=0.0001

        x1=x

        y=0.2

        y1=y

        g=1.0e+37

        dx=1.0

40    u=b(1)

      v=0.0

        do 50 i=1,k

          p=u*x1

          q=v*y1

          pq=(u+v)*(x1+y1)

          u=p-q+b(i+1)

          v=pq-p-q

50    continue

      g1=u*u+v*v

        if(g1.lt.g)goto 105

        if(is.ne.0)goto 80

60    t=t/1.67

      x1=x-t*dx

        y1=y-t*dy

        if(k.ge.50)then

         p=sqrt(x1*x1+y1*y1)

       q=exp(85.0/k)

        if(p.ge.q)goto 60

      end if

        if(t.ge.1.0e-03)goto 40

        if(g.le.1.0e-18)goto 90

65    is=1

      dd=sqrt(dx*dx+dy*dy)

        if(dd.gt.1.0)dd=1.0

        dc=6.28/(k+4.5)

70    c=0.0

80    c=c+dc

      dx=dd*cos(c)

        dy=dd*sin(c)

        x1=x+dx

        y1=y+dy

        if(c.le.6.29)goto 40

        dd=dd/1.67

        if(dd.gt.1.0e-07)goto 70

90    if(abs(y).le.1.0e-06)then

       p=-x

         y=0.0

         q=0.0

        else

         p=-2.0*x

         q=x*x+y*y

         xr(k)=x*w

       xi(k)=-y*w

         k=k-1

        end if

        do 100 i=1,k

         b(i+1)=b(i+1)-b(i)*p

         b(i+2)=b(i+2)-b(i)*p   

 100  continue

      xr(k)=x*w

        xi(k)=y*w

        k=k-1

        if(k.eq.1)then

          xr(k)=-b(2)*w/b(1)

          xi(k)=0.0

          return

        end if

      goto 20

105   g=g1

      x=x1

        y=y1

        is=0

        if(g.le.1.0e-22)goto 90

      u1=k*b(1)

        v1=0.0

        do 110 i=2,k

          p=u1*k

          q=v1*y

          pq=(u1+v1)*(x+y)

          u1=p-q+(k-i+1)*b(i)

          v1=pq-p-q

110   continue

      p=u1*u1+v1*v1

        if(p.le.1.0e-20)goto 65

        dx=(u*u1+v*v1)/p

        dy=(u1*v-v1*u)/p

      t=1.0+4.0/k

        goto 60

        end       

        dimension a(7),xr(6),xi(6),b(7)

        data a/1.0,-5.0,3.0,1.0,-7.0,7.0,-20.0/

        n=7

        m=6

        call dsrrt(a,xr,xi,n,m,l,b)

        if(l.ne.0)then

        do 400 i=1,m

400   continue  

      write(*,500)i,xr(i),xi(i)

      end if

500    format(1x,'x(',i2,1x,')=',e13.6,2x,'j',2x,e13.6)

      end

[ Last edited by 余澤成 on 2010-9-9 at 18:01 ]

回復(fù)此樓

» 猜你喜歡

本人考085602 化學(xué)工程專碩已經(jīng)有23人回復(fù)
求調(diào)劑院校信息已經(jīng)有4人回復(fù)
085600材料與化工306 已經(jīng)有4人回復(fù)
286求調(diào)劑已經(jīng)有10人回復(fù)
328求調(diào)劑，英語六級(jí)551，有科研經(jīng)歷已經(jīng)有9人回復(fù)
一志愿北京化工大學(xué)070300 學(xué)碩336求調(diào)劑已經(jīng)有4人回復(fù)
286分人工智能專業(yè)請(qǐng)求調(diào)劑愿意跨考！已經(jīng)有8人回復(fù)
資源與環(huán)境調(diào)劑申請(qǐng)(333分) 已經(jīng)有5人回復(fù)
280求調(diào)劑已經(jīng)有12人回復(fù)
269專碩求調(diào)劑已經(jīng)有5人回復(fù)

» 本主題相關(guān)價(jià)值貼推薦，對(duì)您同樣有幫助:

fortran程序open文件時(shí)的簡(jiǎn)單問題已經(jīng)有3人回復(fù)
FORTRAN新手求助主程序循環(huán)問題已經(jīng)有10人回復(fù)
Fortran 編譯問題已經(jīng)有9人回復(fù)
我編的Simpson積分法fortran程序給不出結(jié)果，大俠們看看哪里出了問題？已經(jīng)有4人回復(fù)
FORTRAN 基礎(chǔ)知識(shí)講解已經(jīng)有21人回復(fù)
請(qǐng)教一個(gè)fortran小程序編譯出錯(cuò)的問題，謝謝已經(jīng)有9人回復(fù)
Fortran的格式化輸入輸出問題已經(jīng)有14人回復(fù)
有關(guān)fortran的一次而問題，希望大家能幫幫忙，謝謝已經(jīng)有4人回復(fù)
寫了一個(gè)fortran90的小程序，編譯通不過，請(qǐng)大俠幫忙已經(jīng)有59人回復(fù)
【求助】又有新的問題請(qǐng)大家?guī)兔?/a> 已經(jīng)有19人回復(fù)

【交流】Fortran語言答疑專帖已經(jīng)有157人回復(fù)

好好學(xué)習(xí)，天天向上。

1樓 2010-08-07 22:52:12

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

zyj8119

木蟲 (著名寫手)

應(yīng)助: 65 (初中生)
貴賓: 0.003
金幣: 915.1
散金: 1440
紅花: 35
帖子: 2936
在線: 1329.4小時(shí)
蟲號(hào): 664177
注冊(cè): 2008-11-29
性別: GG
專業(yè): 理論和計(jì)算化學(xué)

老是數(shù)組越界，不知道為什么？

贊一下

回復(fù)此樓

好好學(xué)習(xí)，天天向上。

2樓2010-08-07 22:52:38

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

snoopyzhao

至尊木蟲 (職業(yè)作家)

程序強(qiáng)帖: 16
應(yīng)助: 157 (高中生)
貴賓: 0.02
金幣: 18844.7
紅花: 29
帖子: 3803
在線: 1422.4小時(shí)
蟲號(hào): 183750
注冊(cè): 2006-02-13
專業(yè): 污染生態(tài)化學(xué)

★
resonant(金幣+1):感謝提供：-） 2010-08-07 23:55:03
zyj8119(金幣+4):在尋找，你添加的這個(gè)程序好像還是數(shù)組越界。 2010-08-08 16:44:09

問題在于你的數(shù)組 b 下標(biāo)變成負(fù)值了。你跑一下下面的代碼看看，加了幾個(gè) write，可以看出問題出在哪里，呵呵……至于怎么改，我不了解你的算法，沒有辦法……

CODE:

        subroutine dsrrt(a,xr,xi,n,m,l,b)

        dimension a(n),xr(m),xi(m),b(n)

        if(abs(a(1))+1.0.eq.1.0)then

       l=0

         write(*,5)

         return

        end if

5     format(1x,'err')

      l=1

      k=m

      is=0

        w=1.0

        do 10 i=1,n

10    b(i)=a(i)/a(1)

20    pp=abs(b(k+1))

      write(*,*) 'i am at line 16', k, pp

      if(pp.lt.1.0e-12)then

        xr(k)=0.0

        xi(k)=0.0

        k=k-1

        write(*,*) 'i am at line 21', k

        if(k.eq.1)then

         xr(k)=-b(2)*w/b(1)

         xi(k)=0.0

         return

        end if

        goto 20

        end if

        q=pp**(1.0/k)

        p=q

        w=w*p

        do 30 i=1,k

        b(i+1)=b(i+1)/q

        q=q*p

30    continue

      x=0.0001

        x1=x

        y=0.2

        y1=y

        g=1.0e+37

        dx=1.0

40    u=b(1)

      v=0.0

        do 50 i=1,k

          p=u*x1

          q=v*y1

          pq=(u+v)*(x1+y1)

          u=p-q+b(i+1)

          v=pq-p-q

50    continue

      g1=u*u+v*v

        if(g1.lt.g)goto 105

        if(is.ne.0)goto 80

60    t=t/1.67

      x1=x-t*dx

        y1=y-t*dy

        if(k.ge.50)then

         p=sqrt(x1*x1+y1*y1)

       q=exp(85.0/k)

        if(p.ge.q)goto 60

      end if

        if(t.ge.1.0e-03)goto 40

        if(g.le.1.0e-18)goto 90

65    is=1

      dd=sqrt(dx*dx+dy*dy)

        if(dd.gt.1.0)dd=1.0

        dc=6.28/(k+4.5)

70    c=0.0

80    c=c+dc

      dx=dd*cos(c)

        dy=dd*sin(c)

        x1=x+dx

        y1=y+dy

        if(c.le.6.29)goto 40

        dd=dd/1.67

        if(dd.gt.1.0e-07)goto 70

90    if(abs(y).le.1.0e-06)then

       p=-x

         y=0.0

         q=0.0

        else

         p=-2.0*x

         q=x*x+y*y

         xr(k)=x*w

       xi(k)=-y*w

         k=k-1

        end if

        write(*,*) 'i am at line 88', k

        do 100 i=1,k

         b(i+1)=b(i+1)-b(i)*p

         b(i+2)=b(i+2)-b(i)*p   

100   continue

      xr(k)=x*w

        xi(k)=y*w

        k=k-1

        if(k.eq.1)then

          xr(k)=-b(2)*w/b(1)

          xi(k)=0.0

          return

        end if

        write(*,*) 'i am at line 101', k

      goto 20

105   g=g1

      x=x1

        y=y1

        is=0

        if(g.le.1.0e-22)goto 90

      u1=k*b(1)

        v1=0.0

        do 110 i=2,k

          p=u1*k

          q=v1*y

          pq=(u1+v1)*(x+y)

          u1=p-q+(k-i+1)*b(i)

          v1=pq-p-q

110   continue

      p=u1*u1+v1*v1

        if(p.le.1.0e-20)goto 65

        dx=(u*u1+v*v1)/p

        dy=(u1*v-v1*u)/p

      t=1.0+4.0/k

        goto 60

        end      

        program test 

        dimension a(7),xr(6),xi(6),b(7)

        data a/1.0,-5.0,3.0,1.0,-7.0,7.0,-20.0/

        n=7

        m=6

        call dsrrt(a,xr,xi,n,m,l,b)

        if(l.ne.0)then

        do 400 i=1,m

400   continue  

      write(*,500)i,xr(i),xi(i)

      end if

500    format(1x,'x(',i2,1x,')=',e13.6,2x,'j',2x,e13.6)

      end

贊一下(1人)

回復(fù)此樓

3樓2010-08-07 23:31:38

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

weibangxin

★
nono2009(金幣-1):專業(yè)區(qū)請(qǐng)勿純表。謝謝。 2010-08-09 06:43:37
zyj8119(金幣+1): 2010-09-09 08:29:21

4樓2010-08-08 00:10:07

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

snoopyzhao

至尊木蟲 (職業(yè)作家)

程序強(qiáng)帖: 16
應(yīng)助: 157 (高中生)
貴賓: 0.02
金幣: 18844.7
紅花: 29
帖子: 3803
在線: 1422.4小時(shí)
蟲號(hào): 183750
注冊(cè): 2006-02-13
專業(yè): 污染生態(tài)化學(xué)

我說得很清楚啊，我只能找到哪個(gè)地方有問題，但我沒有能力解決問題，呵呵……

贊一下

回復(fù)此樓

5樓2010-08-08 21:25:01

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

snoopyzhao

至尊木蟲 (職業(yè)作家)

程序強(qiáng)帖: 16
應(yīng)助: 157 (高中生)
貴賓: 0.02
金幣: 18844.7
紅花: 29
帖子: 3803
在線: 1422.4小時(shí)
蟲號(hào): 183750
注冊(cè): 2006-02-13
專業(yè): 污染生態(tài)化學(xué)

★
resonant(金幣+1):en,goto太多的程序很煩人奈...哈哈 2010-08-08 23:31:08

試著跟了一下，發(fā)現(xiàn)太困難了，被里面的 goto 搞暈了？這個(gè)程序是你寫的嗎？如果是，想辦法重新理一下結(jié)構(gòu)，盡可能少用 goto……

贊一下(1人)

回復(fù)此樓

6樓2010-08-08 22:27:53

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

zyj8119

木蟲 (著名寫手)

應(yīng)助: 65 (初中生)
貴賓: 0.003
金幣: 915.1
散金: 1440
紅花: 35
帖子: 2936
在線: 1329.4小時(shí)
蟲號(hào): 664177
注冊(cè): 2008-11-29
性別: GG
專業(yè): 理論和計(jì)算化學(xué)

此貼結(jié)貼。

回復(fù)此樓

好好學(xué)習(xí)，天天向上。

7樓2010-09-09 08:29:35

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

最具人氣熱帖推薦 [查看全部]	作者	回/看	最后發(fā)表

[考研] 材料學(xué)碩301分求調(diào)劑 +7	Liyouyumairs 2026-03-21	7/350	2026-03-21 22:31 by peike
[考研] 化學(xué)工程321分求調(diào)劑 +18	大米飯！ 2026-03-15	22/1100	2026-03-21 20:20 by HH領(lǐng)袖
[考研] 考研調(diào)劑 +3	呼呼？~+123456 2026-03-21	3/150	2026-03-21 20:04 by 無際的草原
[考研] 考研化學(xué)學(xué)碩調(diào)劑，一志愿985 +5	張vvvv 2026-03-15	7/350	2026-03-21 19:23 by ColorlessPI
[考研] 26考研一志愿中國(guó)石油大學(xué)(華東)305分求調(diào)劑 +6	嘉年新程 2026-03-15	6/300	2026-03-21 17:07 by Dream007008
[考研] 266求調(diào)劑 +3	哇呼哼呼哼 2026-03-20	3/150	2026-03-21 16:46 by barlinike
[考研] 299求調(diào)劑 +4	某某某某位 2026-03-21	4/200	2026-03-21 16:30 by barlinike
[考研] 求調(diào)劑 +3	.m.. 2026-03-21	4/200	2026-03-21 16:25 by barlinike
[考研] 材料學(xué)學(xué)碩080502 337求調(diào)劑-一志愿華中科技大學(xué) +4	順順順mr 2026-03-18	5/250	2026-03-21 10:22 by luoyongfeng
[考研] 一志愿華南師大 070300（化學(xué)）304分求調(diào)劑 +3	0703武芊慧雪304 2026-03-18	3/150	2026-03-21 00:48 by JourneyLucky
[考研] 一志愿西南交大，求調(diào)劑 +5	材化逐夢(mèng)人 2026-03-18	5/250	2026-03-21 00:26 by JourneyLucky
[考研] 一志愿蘇州大學(xué)材料求調(diào)劑，總分315（英一） +5	sbdksD 2026-03-19	5/250	2026-03-20 22:10 by luoyongfeng
[考研] 260求調(diào)劑 +3	朱芷琳 2026-03-20	3/150	2026-03-20 20:35 by 學(xué)員8dgXkO
[考研] 求調(diào)劑 +3	eation27 2026-03-20	3/150	2026-03-20 19:32 by JourneyLucky
[考研] 281求調(diào)劑（0805） +14	煙汐憶海 2026-03-16	25/1250	2026-03-20 15:47 by yuncha
[考研] 材料與化工求調(diào)劑 +7	為學(xué)666 2026-03-16	7/350	2026-03-19 14:48 by 盡舜堯1
[考研] 材料，紡織，生物（0856、0710），化學(xué)招生啦 +3	Eember. 2026-03-17	9/450	2026-03-18 10:28 by Eember.
[考研] 一志愿南京大學(xué)，080500材料科學(xué)與工程，調(diào)劑 +4	Jy? 2026-03-16	4/200	2026-03-17 11:02 by gaoqiong
[考研] 283求調(diào)劑 +3	聽風(fēng)就是雨； 2026-03-16	3/150	2026-03-17 07:41 by 熱情沙漠
[考研] 085601材料工程315分求調(diào)劑 +3	yang_0104 2026-03-15	3/150	2026-03-15 10:58 by peike