[Wien] Question: geometry optimization: constrains

[B. Yanchitsky]

Peter Blaha wrote:
> WIEN2k cannot do this itself.
> The way to proceed is to change min_lapw and add a small extra program.
> min_lapw extracts the forces and creates case.finM
> 
>   grep :FGL $file.scf | tail -$natom | cut -c29-91 >> $file.finM
> 
> Here you have to insert a small program which reads the forces, and 
> projects out the desired force components and writes them back into 
> case.finM
> In case.inM use NEW1, this should then move the atoms only in the desired 
> plane.
> 
> PS: If you have something more "general", I would be interested to add 
> this to the WIEN2k distribution. 
> 

I have a code for linear constraints, i.e for constraint Ax=b where
A is a (m,3*nat) matrix specified by user, m is number of constraints,
x is 3*nat vectors of atomic coordinates, b is a 3*nat constant vector specified by user.
Each row of matrix A is an equation which defines 3*nat-1 subspace, vector which is normal
to this subspace is parallel to the gradient, but the gradient is just normalised row of matrix A.
We need to substract from the forces their components onto this normal vector and so on for
each row,  after some manipulations it is possible to obtain that constrained forces are
given by
Fconstrained = F - (A^T)*A*F  ,
where ^T means transposition, and A has normalised rows.
I simply added this into mixer which writes constrained forces at the end of scf file,
and added the switch for x_lapw to calculate constrained forces and for min_lapw to grep these forces.
If you are interested in this, i'll send you the file.

ps.
it's somewhat unpractical to change/write even small scripts (programs) for each case,
because this is a good source for errors.

Regards,
Bogdan Yanchitsky,

Institue of Magnetism, Kiev, Ukraine