The problem is not what I posted. After some further testing it appears that the application of the body loading should be done in global frame, thus, in the element 'franf2d.f' (which is what you are using based on your data) you should add
the following before the call to 'fbody2d':
call bm2trn (xl,cs,-sn,ndm*nel,ndm,2)
The file then reads
call bm2trn ( r,cs,-sn,nst,ndf,2)
c Set body loading factors and follower forces
call bm2trn (xl,cs,-sn,ndm*nel,ndm,2)
call fbody2d(d,xl,ul, r,s, ndm,ndf,nst, isw)
The first transforms the residual to the global frame, the second the coordinates and finally computes the body forces.
