Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Small-scale simulation problem 2 #312

Closed
AaronRapopot opened this issue Oct 21, 2023 · 1 comment
Closed

Small-scale simulation problem 2 #312

AaronRapopot opened this issue Oct 21, 2023 · 1 comment
Labels
help wanted Extra attention is needed

Comments

@AaronRapopot
Copy link

Hello, I found that when conducting small-scale simulations, setting the radius too small can lead to numerical errors in the position value. The rod parameters I set refer to "Real-time control of an octopus arm - NNES". I print the position of the rod. The shape of the position is (250, 3, 101). There are a total of 250 steps. In the initial dozens of steps, numerical values are still present (although highly unstable), but all subsequent values displays "nan". If I set the radius from 0.0125 to 0.1, everything is ok.
`
'''Setup Simulation'''
class BaseSimulator(BaseSystemCollection, Constraints, Connections, Forcing, CallBacks, Damping):
pass

'''Create Rod'''
Thermal_drawing_robot = BaseSimulator()

n_elem = 100
L0 = 0.2
radius = 0.012
damp_coefficient = 0.1
rod1 = CosseratRod.straight_rod(
n_elements=n_elem,
start=np.zeros((3,)),
direction=np.array([0.0, 0.0, 1.0]),
normal=np.array([0.0, 1.0, 0.0]),
base_length=L0,
base_radius=radius,
density=700,
youngs_modulus=1e4,
)
Thermal_drawing_robot.append(rod1)

'''Define Boundary Conditions'''
Thermal_drawing_robot.constrain(rod1).using(
OneEndFixedRod,
constrained_position_idx=(0,),
constrained_director_idx=(0,)
)

'''Define torque'''
torque = 1.0
torque_direction = np.array([0.0, -1.0, 0.0])
Thermal_drawing_robot.add_forcing_to(rod1).using(
UniformTorques, torque, direction=torque_direction
)

'''Apply damping'''
time_step = 1.0e-5
nu = damp_coefficient
Thermal_drawing_robot.dampen(rod1).using(
AnalyticalLinearDamper,
damping_constant = nu,
time_step = time_step,
)

'''Add callback'''
class RodCallBack(CallBackBaseClass):
def init(self, step_skip: int, callback_params: dict):
CallBackBaseClass.init(self)
self.every = step_skip
self.callback_params = callback_params

def make_callback(self, system, time, current_step: int):
    if current_step % self.every == 0:
        self.callback_params["time"].append(time)
        self.callback_params["radius"].append(system.radius.copy())
        self.callback_params["position"].append(system.position_collection.copy())
        self.callback_params["director"].append(system.director_collection.copy())
        self.callback_params["velocity"].append(system.velocity_collection.copy())
        return

callback_data_rod1 = defaultdict(list)
Thermal_drawing_robot.collect_diagnostics(rod1).using(
RodCallBack, step_skip=1000, callback_params=callback_data_rod1)

'''Finalize simulator'''
Thermal_drawing_robot.finalize()

'''Set timestepper'''
timestepper = PositionVerlet()
final_time = 2.5 # seconds
total_steps = int(final_time / time_step)
integrate(timestepper, Thermal_drawing_robot, final_time, total_steps)

print(callback_data_rod1['position'])
position_shape = np.array(callback_data_rod1["position"]).shape
print(position_shape)
And the result (just show one step, others are same):
array([[ 0., nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan],
[ 0., nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan],
[ 0., nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
nan, nan, nan, nan, nan, nan, nan, nan, nan, nan]])
`

@AaronRapopot AaronRapopot added the help wanted Extra attention is needed label Oct 21, 2023
@AaronRapopot
Copy link
Author

I solved it. Set youngs modulus from 1e4 to 1e7, it works.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
help wanted Extra attention is needed
Projects
None yet
Development

No branches or pull requests

2 participants