Dear Professor Taylor,
Thank you very much for your response. That solved the problem. Just one comment: apparently somewhere between the main feap and nurbfeap, hsize array is restarted because when I write out the hsize array in plot routines in nurbfeap (e.g. pmsh_lin.f ), it returns a zero array. Fortunately, calling phsize again within the nurbfeap module recalculates the hsize and then the plot routines work perfectly (currently, I've added the call to the beginning of pform.f, since I think it is one of the first routines which feap calls).
BTW, I would like to thank you once again for your patience during this debugging process.