I checked some of the problem cases where I struggled in the past with the “RF” output. Some of them dealt with gravity loads (e.g. CCX gravity problem? ). For the mentioned gravity problem example the output for RF and RR is equal because there is no point load, only a distributed load through *DLOAD. Is that what you want to achieve through your development as you mentioned in your initial post in this thread to adress the point loads only?
Here the box gravity model with some additional output like mass and RR:
*node
1, 0., 0., 0.
2, 10., 0., 0.
3, 10., 10., 0.
4, 0., 10., 0.
**
5, 0., 0., 10.
6, 10., 0., 10.
7, 10., 10., 10.
8, 0., 10., 10.
**
9, 5., 0., 0.
10, 10., 5., 0.
11, 5., 10., 0.
12, 0., 5., 0.
**
13, 5., 0., 10.
14, 10., 5., 10.
15, 5., 10., 10.
16, 0., 5., 10.
**
17, 0., 0., 5.
18, 10., 0., 5.
19, 10., 10., 5.
20, 0., 10., 5.
**
*nset, nset=bottom
1, 2, 5, 6, 9, 13, 17, 18
**
*element, type=c3d20, elset=cube
1, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15
16, 17, 18, 19, 20
**
*material, name=steel
*density
7.8E-09
*elastic
210000, 0.28
**
*solid section, elset=cube, material=steel
**
** +++++++++++++++++++++++++++++++++++++++
*step
*static
1.,1.,1.E-05,1.
**
*boundary, op=new
bottom, 1, 3, 0
**
*Dload
cube, Grav, -9800., 0.0, 1.0, 0.0
**
*node print, nset=bottom, totals=Yes
rf, RR
**
*node file
rf, u, RR
*el file
s, e
*EL PRINT, TOTALS=YES, ELSET=CUBE
EVOL, EMAS
*end step