fix addforce and atom variable

Hi all,

I’m trying to add a drag force on each particles. The drag is caused by a mean “fluid” velocity that varies linearly with height.
I did :

Hi all,

I'm trying to add a drag force on each particles. The drag is caused by a
mean "fluid" velocity that varies linearly with height.
I did :

_____________________
# Constants
variable stress equal 7
variable Cd equal 0.41
variable z0 equal 0.04

# Get atom radius
compute d all property/atom radius

#Velocity profile
variable Ux atom (v_z-z0)*stress

#Drag force in x direction on particle
variable fxdrag atom 1000/2*(2*c_d)^2*Cd*(v_Ux-v_vx)*(v_Ux-v_vx)

note: you use Cd here and not v_Cd

fix kick all addforce v_fxdrag 0 0
_____________________

However it does not work properly. When I try for instance to print the
average drag force on particle with :
_______________________________
compute meandrag all reduce ave v_fxdrag
variable md equal c_meandrag
fix extra all print 100 "\{time\} {md}" file ./post/meandrag.txt
title "%" screen yes
________________________________

I get an error at the print level : " ERROR on proc 0: Substitution for
illegal variable".
I tried even more simple, to print the mean radius of particle with :

_______________________________
compute meanrad all reduce ave c_d
variable md equal c_meanrad
fix extra all print 100 "\{time\} {md}" file ./post/meandrag.txt
title "%" screen yes
________________________________

but I get also "ERROR on proc 0: Substitution for illegal variable".

I guess I'm not defining correctly my variables... Some suggestions ?

the way to go about debugging any of these is
to remove *all* variable related lines and build them
line by line and step by step. the error you are seeing
is usually due to using a string where you forgot the
v_ or c_ prefix, so that lammps will try to interpret it
as a "thermo" property, but not on the top level, but
inside a variable definition, causing the variable
definition to fail.

axel.

Thks alex for this quick reply.

Cd is a constant in my case (drag coefficient) while c_d should give the vector of particles radius.

For the debugging, the funny thing is that without the "print" lines but with all variable definition, there is no error. But I'd like to be sure of the force values applied, so I wanted to print it.

An other question: In the addforce fix, I give a vector for the fx force (corresponding to the forces applied to each particles v_fxdrag) and not a constant, as often it is the case in examples. Is that suppose to work fine too ?

Thks a lot,
j

Thks alex for this quick reply.

Cd is a constant in my case (drag coefficient) while c_d should give the
vector of particles radius.

constant or not, you have to refer to it either via v_Cd or ${Cd}

For the debugging, the funny thing is that without the "print" lines but
with all variable definition, there is no error. But I'd like to be sure of
the force values applied, so I wanted to print it.

fix addforce may not be able to propagate error messages.
using variables in it was added after the fact, so there is
always a chance of unexpected side effects.

getting the print to work is a good strategy.

An other question: In the addforce fix, I give a vector for the fx force
(corresponding to the forces applied to each particles v_fxdrag) and not a
constant, as often it is the case in examples. Is that suppose to work fine
too ?

yes, atom style variables are supported.
it should say so in the documentation.

axel.

Ok for Cd.

Let's take a simpler case. If I do:

variable Cd equal 0.41
fix kick all addforce ${Cd} 0 0

It adds correctly the force.

But if I do:

variable Cd *atom* 0.41
fix kick all addforce v_Cd 0 0

The force is no longer transmitted, even if it should be the same (vector of Cd for each atom, or constant Cd)...

I'm working with LIGGGHTS 1.5 based on lammps-10Mar10, could it be a problem ?

j

Ok for Cd.

Let's take a simpler case. If I do:

variable Cd equal 0.41
fix kick all addforce ${Cd} 0 0

It adds correctly the force.

But if I do:

variable Cd atom 0.41
fix kick all addforce v_Cd 0 0

The force is no longer transmitted, even if it should be the same (vector of
Cd for each atom, or constant Cd)...

I'm working with LIGGGHTS 1.5 based on lammps-10Mar10, could it be a problem
?

that *is* the problem. you have to refer to the
documentation of your version. the online
version always refers to the current version.
LAMMPS has changed a lot since 2010.
using per atom variables for fix addforce/setforce
is one of the added features.

in general, if you don't specify which version
you use, it is always implied that you use
the latest version of LAMMPS.

axel.