求助，子程序中一个自定义函数的问题

Harry

[Fortran] 纯文本查看 复制代码

查看代码

打印?

001	subroutine Front

002	dimension fixed(160),equat(60),vecrv(160),gload(60),gstif(1830),estif(16,16),&

003	iffix(160),nacva(60),locel(16),ndest(16)

004	common/Contro/npoin,nelem,nnode,ndofn,ndime,nstre,ntype,ngaus,nprop,nmats,nvfix,nevab,&

005	icase,ncase,itemp,iprob,nprob

006	common/Lgdata/coord(80,2),props(10,5),presc(40,2),asdis(160),eload(25,16),strin(3,225),&

007	nofix(40),ifpre(40,2),lnods(25,8),matno(25)

008	common/Work/elcod(2,8),shape(8),deriv(2,8),dmatx(3,3),cartd(2,8),dbmat(3,16),bmatx(3,16),&

009	smatx(3,16,9),posgp(3),weigp(3),gpcod(2,9),neror(24)

010	!integer nfunc

011

012	call Openfile

013

014

015	nfunc(i,j)=(j*j-j)/2+i            !*********************************问题在此行

016

017

018

mfron=60

019	mstif=1830

020

!

021	!  interpret fixity data in vector form

022

!

023	ntotv=npoin*ndofn

024	do 100 itotv=1,ntotv

025	iffix(itotv)=0

026	100 fixed(itotv)=0.0

027	do 110 idofn=1,ndofn

028	ngash=nloca+idofn

029	iffix(ngash)=ifpre(ivfix,idofn)

030	110 fixed(ngash)=presc(ivfix,jdofn)

031

!

032	!  change the sign of the last appearance

033	!  of each node

034

!

035	do 140 ipoin=1,npoin

036

klast=0

037	do 130 ielem=1,nelem

038	do 120 inode=1,nnode

039	if(lnods(ielem,inode).ne.ipoin) goto 120

040	klast=ielem

041	nlast=inode

042	120 continue

043	130 continue

044	if(klast.ne.0) lnods(klast,nlast)=-ipoin

045	140 continue

046

!

047	!  start by initializing everything that

048	!  matters to zero

049

!

050	do 150 istif=1,mstif

051	150 gstif(istif)=0.0

052	do 160 ifron=1,mfron

053	gload(ifron)=0.0

054	equat(ifron)=0.0

055	vecrv(ifron)=0.0

056	160 nacva(ifron)=0

057

!

058	!  and prepare for disc reading and writing

059	!  operations

060

!

061

rewind 2

062

rewind 4

063

rewind 6

064

rewind 8

065

!

066	!  enter main element assembly-reduction loop

067

!

068

nfron=0

069

kelva=0

070	do 380 ielem=1,nelem

071

kevab=0

072	read(2,*) estif

073	do 170 inode=1,nnode

074	do 170 idofn=1,ndofn

075	nposi=(inode-1)*ndofn+idofn

076	locno=lnods(ielem,inode)

077	if(locno.gt.0) locel(nposi)=(locno-1)*ndofn-idofn

078	if(locno.lt.0) locel(nposi)=(locno+1)*ndofn-idofn

079	170 continue

080

!

081	!  start by looking for existing destinations

082

!

083	do 210 ievab=1,nevab

084	nikno=iabs(locel(ievab))

085

kexis=0

086	do 180 ifron=1,nfron

087	if(nikno.ne.nacva(ifron)) goto 180

088	kevab=kevab+1

089

kexis=1

090	ndest(kevab)=ifron

091	180 continue

092	if(kexis.ne.0) goto 210

093

!

094	!  we now seek new empty places for

095	!  destination vector

096

!

097	do 190 ifron=1,mfron

098	if(nacva(ifron).ne.0) goto 190

099	nacva(ifron)=nikno

100	kevab=kevab+1

101	ndest(kevab)=ifron

102	190 continue

103

!

104	!  the new places may demand an increase

105	!  in current frontwidth

106

!

107	200 if(nedst(kevab).gt.nfron) nfron=ndest(kevab)

108	210 continue

109

!

110	!  assemble element loads

111

!

112	do 240 ievab=1,nevab

113	idest=ndest(ievab)

114	gload(idest)=gload(idest)+eload(ielem,ievab)

115

!

116	!  assemble the element stifnesses

117	!  but not in resolution

118

!

119	if(icase.gt.1) goto 230

120	do 220 jevab=1,ievab

121	jdest=ndest(jevab)

122	ngash=nfunc(idest,jdest)

123	ngish=nfunc(jdest,idest)

124	if(jdest.ge.idest) gstif(ngash)=gstif(ngash)+estif(ievab,jevab)

125	if(jdest.lt.idest) gstif(ngish)=gstif(ngish)+estif(ievab,jevab)

126	220 continue

127	230 continue

128	240 continue

129

!

130	!  re-examine each element node, to

131	!  enquire which can be eliminated

132

!

133	do 370 ievab=1,nevab

134	nikno=-locel(ievab)

135	if(nikno.le.0) goto 370

136

!

137	!  find positions of variables ready

138	!  for elimination

139

!

140	do 350 ifron=1,nfron

141	if(nacva(ifron).ne.nikno) goto 350

142

!

143	!  extract the cofficients of the

144	!  new equation for elimination

145

!

146	if(icase.gt.1) goto 260

147	do 250 jfron=1,mfron

148	if(ifron.lt.jfron) nloca=nfunc(ifron,jfron)

149	if(ifron.ge.jfron) nloca=nfunc(jfron,ifron)

150	equat(jfron)=gstif(nloca)

151	250 gstif(nloca)=0.0

152	260 continue

153

!

154	!  and fxtract the corresponding right

155	!  hand sides

156

!

157	eqrhs=gload(ifron)

158	gload(ifron)=0.0

159	kelva=kelva+1

160

!

161	!  write equations to disc or to tape

162

!

163	if(icase.gt.1) goto 270

164	write(4,*) equat,eqrhs,ifron,nikno

165

goto 280

166	270 write(8,*) eqrhs

167	read(4,*) wquat,dummy,idumm,nikno

168	280 continue

169

!

170	!  deal with pivot

171

!

172	pivot=equat(ifron)

173	equat(ifron)=0.0

174

!

175	!  enquire whether present varizable is

176	!  free or prescribed

177

!

178	if(iffix(nikno).eq.0) goto 300

179

!

180	!  deal with a prescribed deflection

181

!

182	do 290 jfron=1,nfron

183	290 gload(jfron)=gload(jfron)-fixed(nikno)*equat(jfron)

184

goto 340

185

!

186	!  eliminate a free variable-deal with

187	!  the right hand side first

188

!

189	300 do 330 jfron=1,nfron

190	gload(jfron)=gload(jfron)-equat(jfron)*eqrhs/pivot

191

!

192	!  now deal with the coefficients in core

193

!

194	if(icase.gt.1) goto 320

195	if(equat(jfron).eq.0.0) goto 330

196	nloca=nfunc(0,jfron)

197	do 310 lfron=1,jfron

198	ngash=lfron+nloca

199	310 gstif(ngash)=gstif(ngash)-equat(jfron)*equat(lfron)/pivot

200	320 continue

201	330 continue

202	340 equat(ifron)=pivot

203

!

204	!  record the new vacant space, and reduce

205	!  frontwidth if possible

206

!

207	nacva(ifron)=0

208

goto 360

209

!

210	!  complete the element loop in the forward

211	!  elimination

212

!

213	350 continue

214	360 if(nacva(nfron).ne.0) goto 370

215	nfron=nfron-1

216	if(nfron.gt.0) goto 360

217	370 continue

218	380 continue

219

!

220	!  enter back-substitution phase, loop

221	!  backwards through variables

222

!

223	do 410 ielva=1,kelva

224

!

225	!  read a new equation

226

!

227	backspace 4

228	read(4,*) equat,eqrhs,ifron,nikno

229	backspace 4

230	if(icase.eq.1) goto 390

231	backspace 8

232	read(8,*) eqrhs

233	backspace 8

234	390 continue

235

!

236	!  prepare to back-substitute from the

237	!  current equation

238

!

239	pivot=equat(ifron)

240	if(iffix(nikno).eq.1) vecrv(ifron)=fixed(nikno)

241	if(iffix(nikno).eq.0) equat(ifron)=0.0

242

!

243	!  back-substitute in the current equation

244

!

245	do 400 jfron=1,mfron

246	400 eqrhs=eqrhs-vecrv(jfron)*equat(jfron)

247

!

248	!  put the final values where they belong

249

!

250	if(iffix(nikno).eq.0) vecrv(ifron)=eqrhs/pivot

251	if(iffix(nikno).eq.1) fixed(nikno)=-eqrhs

252	asdis(nikno)=vecrv(ifron)

253	410 continue

254	write(12,900)

255	900 format(1ho,5x,13hDisplacements)

256	if(ndofn.ne.2) goto 430

257	if(ndime.ne.1) goto 420

258	write(12,905)

259	905 format(1ho,5x,4hnode,6x,5hDisp.,7x,8hRotation)

260

goto 440

261	420 write(12,910)

262	910 format(1ho,5x,4hNode,5x,7hX-disp.,7x,7hY-disp.)

263	430 write(12,915)

264	915 format(1ho,5x,4hNode,5x,7hX-disp.,7x,7hY-disp.)

265	440 continue

266	do 450 ipoin=1,npoin

267	ngash=ipoin*ndofn

268	ngish=ngash-ndofn+1

269	450 write(12,920) ipoin,(asdis(igash),igash=ngish,ngash)

270	920 format(i10,3e14.6)

271	write(12,925)

272	925 format(1ho,5x,9hReactions)

273	if(ndofn.ne.2) goto 470

274	if(ndime.ne.1) goto 460

275	write(12,930)

276	930 format(1ho,5x,4hNode,6x,5hForce,8x,6hMomnet)

277

goto 480

278	460 write(12,935)

279	935 format(1ho,5x,4hNode,5x,7hX-force,7x,7hY-froce)

280

goto 480

281	470 write(12,940)

282	940 format(1ho ,5x,4hNode,6x,5hForce,6x,9hXz-moment,5x,9hYz-moment)

283	480 continue

284	do 510 ipoin=1,npoin

285	nloca=(ipoin-1)*ndofn

286	do 490 idofn=1,ndofn

287	ngush=nloca+idofn

288	if(iffix(ngush).gt.0) goto 500

289	490 continue

290

goto 510

291	500 ngash=nloca+ndofn

292	ngish=nloca+1

293	write(12,945) ipoin,(fixed(igash),igash=ngish,ngash)

294	510 continue

295	945 format(i10,3e14.6)

296

!

297	!  post front-reset all element connection numbers to positive

298	!  values for subsequent use in stress calculation

299

!

300	do 520 ielem=1,nelem

301	do 520 inode=1,nnode

302	520 lnods(ielem,inode)=iabs(lnods(ielem,inode))

303

return

304

305

end

vvt

您有何问题？请给出文字描述。

fcode

nfunc(i,j)=(j*j-j)/2+i            !*********************************问题在此行

移动到 call Openfile 前面

Harry

移了，但还是有问题
Error        1         error LNK2019: unresolved external symbol _NEDST referenced in function _FRONT        s-front.obj       
Error        2         fatal error LNK1120: 1 unresolved externals        Debug\master-plane.exe       

fcode

第 107 行
200 if(nedst(kevab).gt.nfron) nfron=ndest(kevab)
是否应该是
200 if(ndest(kevab).gt.nfron) nfron=ndest(kevab)

Harry

是的，找出错误了，太感谢了！

百事可乐

implicit none 你值得拥有。

Harry

fcode 发表于 2015-12-22 08:46
第 107 行
200 if(nedst(kevab).gt.nfron) nfron=ndest(kevab)
是否应该是

大侠，编译通过了，可是，运行过程中当call Front之前，整形数组nofix(40),ifpre(40,2),lnods(25,8),matno(25)是有内容的，可是进入Front后这些数组全部回到了0值，其他非整形的数组没有问题。不知道是什么原因？？？（openfile中没有用implicit none）

fcode

不知道，运行以后的问题属于动态问题，需要动态分析。（调试等）

而只有代码的一部分，无法进行动态分析。