<div dir="ltr"><div>Dear Dorothea,</div><div>Thank you for your answer. </div><div> </div><blockquote class="gmail_quote" style="margin: 0;margin-left: 0.8ex;border-left: 1px #ccc solid;padding-left: 1ex;"><div dir="ltr"><div>I am not sure if I understand your question correctly </div></div></blockquote><div> </div><div>On the longer run, we are interested in describing an electrolyte above a metal surface and using a charge imbalance to create a potential drop at the interface very much as one would do in a pure DFT calculation, but then withouth the need to describe the metal quantum mechanically. I know that this is not the original idea behind the image charge method as implemented in cp2k, but I thought it may offer a work around. So to understand how things work, I simply considered an H2+ molecule above a surface. This is similar to the model system used by Siepmann and Sprik in their 1995 publication, just that they used an MM point charge above a metal. In principle, it should be possible to determine V_0 from the condition that the energy should be minimized, and the resulting counter charge should exactly cancel the charge in front of the metal. </div><div> </div><blockquote class="gmail_quote" style="margin: 0;margin-left: 0.8ex;border-left: 1px #ccc solid;padding-left: 1ex;"><div dir="ltr"><div>and what your computational setup is.</div></div></blockquote><div>I paste here an example of the input that I used with MT boundary conditions. The resulting plot for the Hartree-potential (obtained using &V_HARTREE_CUBE and then averaging over x and y to get a plot over z) is shown in the third image in my original post (green line). </div><div>As you can see, I simply adapted the input file from the GUA example on the web, by replacing the GUA in the xyz file by H2, setting the DFT charge to 1 and changing the poisson solver in the DFT and MM part. The potential EXT_POTENTIAL was determined by running 2 calculations with changing EXT_POTENTIAL_1 and EXT_POTENTIAL_2, and then determining EXT_POTENTIAL_3 by requiring that the resulting charge should be -1 and assuming a linear relationship (which for this system is perfectly followed, giving me a total charge on the metal of -1)</div><div> </div><div>-----------------Begin sample input file using the MT poisson solver------------------ </div><div><div style="background-color: rgb(250, 250, 250); border-color: rgb(187, 187, 187); border-style: solid; border-width: 1px; overflow-wrap: break-word;" class="prettyprint"><code class="prettyprint"><div class="subprettyprint">&FORCE_EVAL METHOD QMMM &DFT CHARGE 1 LSD BASIS_SET_FILE_NAME BASIS_MOLOPT POTENTIAL_FILE_NAME GTH_POTENTIALS &POISSON &EWALD EWALD_TYPE ewald ALPHA .44 GMAX 21 &END EWALD POISSON_SOLVER MT PERIODIC XY &END POISSON &MGRID COMMENSURATE CUTOFF 300 NGRIDS 5 &END MGRID &QS METHOD GPW EXTRAPOLATION ASPC EXTRAPOLATION_ORDER 3 &END QS &SCF MAX_SCF 300 SCF_GUESS ATOMIC # SCF_GUESS RESTART EPS_SCF 1.0E-6 &OT PRECONDITIONER FULL_SINGLE_INVERSE MINIMIZER DIIS &END &END SCF &XC &XC_FUNCTIONAL PBE &END XC_FUNCTIONAL &VDW_POTENTIAL DISPERSION_FUNCTIONAL PAIR_POTENTIAL &PAIR_POTENTIAL TYPE DFTD3 CALCULATE_C9_TERM .TRUE. REFERENCE_C9_TERM PARAMETER_FILE_NAME dftd3.dat REFERENCE_FUNCTIONAL PBE R_CUTOFF [angstrom] 16.0 &END PAIR_POTENTIAL &END VDW_POTENTIAL &END XC &PRINT &V_HARTREE_CUBE &END &END &END DFT &MM &FORCEFIELD &CHARGE ATOM Au CHARGE 0 &END CHARGE &CHARGE ATOM H CHARGE 0 &END CHARGE &SPLINE EPS_SPLINE 1.E-5 #EMAX_SPLINE 2.0 &END &NONBONDED &EAM atoms Au Au PARM_FILE_NAME Au.pot &END EAM &LENNARD-JONES atoms Au H EPSILON 0.0 SIGMA 3.166 RCUT 15 &END LENNARD-JONES &LENNARD-JONES atoms H H EPSILON 0.0 SIGMA 3.166 RCUT 15 &END LENNARD-JONES &END &END FORCEFIELD &POISSON &EWALD EWALD_TYPE ewald ALPHA .44 GMAX 21 &END EWALD POISSON_SOLVER MT PERIODIC XY &END POISSON &END MM &QMMM CENTER NEVER &CELL ABC [angstrom] 34.6055 29.9693 60.0 PERIODIC XY &END CELL &QM_KIND H MM_INDEX 577..578 &END QM_KIND &IMAGE_CHARGE EXT_POTENTIAL 1.24534 MM_ATOM_LIST 1..576 WIDTH 3.5 &END IMAGE_CHARGE &PRINT &IMAGE_CHARGE_INFO &END &END &END QMMM &SUBSYS &CELL ABC [angstrom] 34.6055 29.9693 60.0 PERIODIC XY &END CELL &TOPOLOGY COORD_FILE_NAME Au_gua_image_dampFunc-opt.xyz COORDINATE xyz &END &KIND H BASIS_SET DZVP-MOLOPT-SR-GTH POTENTIAL GTH-PBE-q1 &END KIND &END SUBSYS &END FORCE_EVAL &GLOBAL PROJECT Au_gua_image_dampFunc RUN_TYPE energy &END GLOBAL</div></code></div> </div><div>-----------------End sample input file using the MT poisson solver------------------</div><div> </div><div>>If you have a charged system in the QM part, then you need a solver like MT, like in a standard QM calculation.</div><div>Yes, I would have thought, that this is what I am doing, but the resulting plot for the hartree potential looks strange (shielding charges in the metal, leading to a potential gradient in the metal instead of a constant) and I cannot seem to find an energy minimum, when the potential is chosen such, that the image charge is exatly shielded, which was one of the test systems in the Siepmann Sprik paper (https://aip.scitation.org/doi/abstract/10.1063/1.469429).</div><div> </div><div>I guess the latter is my main issue, as we will depend heavily on the energies. (Hence I need to understand what they contain)</div><div> </div><div>Thank you for your help,</div><div>Katharina</div><div> </div><div> </div> </div>