Dear all,
I'm having what seemed like a trivial issue when applying fix addforce to individual groups of atoms, but in spite of trying for half a day to solve it I cannot seem to find a solution, so I was hoping someone else might be able to help!
What I'm trying to do is to use fix addforce to add a "propulsion force" to each particle, which should at each instant be parallel to its dipole vector. However, I also want the magnitude of the force to vary based on the particle's location; this is done as usual by defining regions and then a dynamic group assignment updated every 100 timesteps. The force is then calculated and added through the following syntax:
compute myMu all property/atom mux muy muz
variable force_x atom c_myMu[1] # Force is parallel to dipole vector (which has magnitude 1)
variable force_y atom c_myMu[2]
variable force_z atom c_myMu[3]
variable f0_x atom v_f0*v_force_x # f0, f1, f2... are constants defined earlier
variable f0_y atom v_f0*v_force_y
variable f0_z atom v_f0*v_force_z
fix propForce0 r0 addforce v_f0_x v_f0_y v_f0_z # Groups are named r0, r1, r2...
variable f1_x atom v_f1*v_force_x
variable f1_y atom v_f1*v_force_y
variable f1_z atom v_f1*v_force_z
fix propForce1 r1 addforce v_f1_x v_f1_y v_f1_z
.
So far, so good! The problem is that I have *17* different regions where I want different force magnitudes, and defining atom-style variables f1_x, f2_x, etc., for all these 17 groups is a very inefficient way of doing things, since it means I'm calculating 17 different forces for each atom on every timestep, when I only actually need one. What I would need is a single atom style variable where the individual element gives the magnitude of the propulsion force for that particular atom, based on the current group assignment.
One alternative solution would be to perform the multiplication directly inside the addforce command:
fix propForce0 v0 addforce \(v\_f0\*v\_force\_x\) (v_f0*v_force_y) $(v_f0*v_force_z)
This however generates an error message that the force variable is of style "equal" rather than "atom", which is what is required. I have also tried to find a solution based on fix property/atom, but it doesn't seem to be viable either.
Using the dirty solution above is hogging down the performance and scalability to the limit of being unusable, which is not surprising. Does anyone know a way around this without hacking the code? And, if not, what do you think would be the easiest way of implementing it into the code?
Any input is welcome, thank you!
Best regards,
Joakim Stenhammar