rotation in local mode

My thanks to everyone who has replied to this question. However, I am going to rewrite this as it is obvious
that I have not made the question clear enough - my apologies.

Let's say I rez 4 cubes in front of me (the actual number is immaterial and they are not linked) and I want to send
them all a message to rotate 25 degrees in the world x-axis, but at the same time keeping their relative positions.

I can do this in a similar fashion to the answer give by Void Singer.

1) I find the centre of the 4 objects and calculate each object's offset from that centre.

2) I add the 25 degrees rotation to each object (their rotation + 25 degrees)

3) Their new position will have changed by their offset multiplied by the added rotation.

4) I add this to their original position.

So far so good - they will rotate 25 degrees and keep their relative positions. But of course, when I originally
rezzed them they had zero rotation and lie firmly in the world axes.

Now let's say I turn 45 degrees and move the cubes in front of me as before (ie. rotating them 45 degrees). If
I now send them a message they will rotate and keep their relative postions but their rotation will be a mixture
of x-axis and y-axis because they are now longer aligned with the world axes.

To remedy this I can reverse the rotation in 2) above: instead of rotation + 25 degrees I can put 25 degrees +
rotation. They will rotate on their local axes (which is what I want) but their relative positions will be
compromised.

The question is: how do you calulate the new positions so that they remain relative to their original positions?

I think you are saying they are linked together, right?

If so, then the three children will have a new position and rotation reletaive to the parent. You will need to calculate those.

There is a tool called Puppeteer that might help. You set up the object, record a snapshot, then move all the parts, another snapshot....etc.

There are scripts that record everything, and then play it all back.

This is pseudo-code, but you should get the idea:

vector CENTER_OF_ROTATION = ...;

// Axial direction AND radians/sec as magnitude, or split to optimize as you will
vector OMEGA = ...;

float initialTime = llGetTime();
vector initialDisplacement = llGetPos()-CENTER_OF_ROTATION;
rotation initialRot = llGetRot();

...

float elapsedTime = llGetTime()-initialTime;
float theta = llVecMag(OMEGA)*elapsedTime;
vector axis = llVecNorm(OMEGA);
rotation deltaRot = llAxisAngle2Rot(axis, theta);

vector pos = CENTER_OF_ROTATION+initialDisplacement*deltaRot;
rotation rot = initialRot*deltaRot;

llSetPrimitiveParams([ PRIM_POSITION, pos, PRIM_ROTATION, rot]);
// Or some equivalent llMoveToTarget() and llRotLookAt() combination
// if physical, but this takes thinking a little more dynamically unless you do it very often.

4 prims not touching, but part of the same object?

rotating the whole object, around it's root prim is straightforward, just rotate the root, the child prims will match it's orientation (rotation, distance and direction).

rotating the child prims in place from a script in that prim
llSetLocalRot( rot );
//-- below are workarounds because of a bug,
llSetRot( rot / llGetRootRotation() );
llSetPrimitiveParams( [PRIM_ROTATION, rot / llGetRootRotation()] );

from a script in another prim
//-- same workaround
llSetLinkPrimitiveParams( <insert link number>, [PRIM_ROTATION, rot / llGetRootRotation()] );

rotating the root without rotating the rest of the set
PITA, you must move the position of all other prims to represent the opposite change in rot as the root, as well as change their rotation that amount, and then rotate the root prim... multiple steps

Emphasis added:

just covering the bases, since OP description isn't 100% clear (object implies linking), figured they might have meant 'not touching'.

if they are truly unlinked, and OP wants to keep their positions the same, but not their individual rotations, they'd be rotated around their individual local axis by
llSetRot( <amount to rotate> * llGetRot() );

or their individual global axis with
llSetRot( llGetRot() * <amount to rotate> );

another possibility is getting the center point of the objects, either as a known offset, or an average of the objects positions, and rotating them around that point similar to what you did only I have a slightly different formula for position and offset
offset = <average all positions>; //-- or average local offsets using 0 for the root, and ad to get pos for a linked object
arc2travel = llEuler2Rot( <degrees to change> * DEG_TO_RAD ); //-- or a rot to add to the current rotation;
pos = llGetPos() + (offset - offset * arc2travel) * llGetRot();
rot = arc2travel * llGetRot();