### Godot version

4.0.2

### System information

Windows 10 960m

### Issue description

When using the skeltonIK3d on a simple bone chain, the roll axis of bones are not working correctly. As can be seen in this clip comparing 3.5 to 4.0.2, when moving the ik target on the x axis, the bone roll of the shin bone twists in a very strange way.

When the ik target is just moving up and down on the y axis 4.0.2 looks the same as 3.5. I have tried many ways of setting the roll axis in blender before exporting, but this seems to make no difference in godot.

https://github.com/godotengine/godot/assets/7741669/524f9450-8582-4336-96a2-397ffe915d2d

@lyuma I also tried this with your recent fix in master and it does not fix this issue, unfortunately. I did try the scorpion demo scene and I can see that it resolved that. Simple 2 or 3 bone chains seem to exhibit this odd roll issue...maybe the long chain on the scorpion tail is covering up the issue?

On a side note, I am getting tons of crashes trying to test out simple ik with the standard tool script. Especially just clicking in the skeletonik3d node. 

### Steps to reproduce

- import a simple 2 or 3 bone skeleton such as a leg
- add a skeltonik3d node
- fill in the node parameters and add a marker3d as target
- add a tool script to start off the ik solve
- move the ik target on the x axis and observe the twisting

### Minimal reproduction project

This has the 3.5, 4.0.2 file and the blend file

[ik_roll_bug.zip](https://github.com/godotengine/godot/files/11521836/ik_roll_bug.zip)

}
		r_axis = Vector3(x, y, z);
		r_angle = Math_PI;
		return;
	}
	
	double s = Math::sqrt((rows[2][1] - rows[1][2]) * (rows[2][1] - rows[1][2]) + (rows[0][2] - rows[2][0]) * (rows[0][2] - rows[2][0]) + (rows[1][0] - rows[0][1]) * (rows[1][0] - rows[0][1])); 
	if (Math::abs(s) < CMP_EPSILON) {
		
		s = 1;
	}
	x = (rows[2][1] - rows[1][2]) / s;
	y = (rows[0][2] - rows[2][0]) / s;
	z = (rows[1][0] - rows[0][1]) / s;
	r_axis = Vector3(x, y, z);
	
	r_angle = Math::acos((rows[0][0] + rows[1][1] + rows[2][2] - 1) / 2);
}
void Basis::set_quaternion(const Quaternion &p_quaternion) {
	real_t d = p_quaternion.length_squared();
	real_t s = 2.0f / d;
	real_t xs = p_quaternion.x * s, ys = p_quaternion.y * s, zs = p_quaternion.z * s;
	real_t wx = p_quaternion.w * xs, wy = p_quaternion.w * ys, wz = p_quaternion.w * zs;
	real_t xx = p_quaternion.x * xs, xy = p_quaternion.x * ys, xz = p_quaternion.x * zs;
	real_t yy = p_quaternion.y * ys, yz = p_quaternion.y * zs, zz = p_quaternion.z * zs;
	set(1.0f - (yy + zz), xy - wz, xz + wy,
			xy + wz, 1.0f - (xx + zz), yz - wx,
			xz - wy, yz + wx, 1.0f - (xx + yy));
}
void Basis::set_axis_angle(const Vector3 &p_axis, real_t p_angle) {