Search in sources :

Example 66 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class ParticleMeshEwaldQI method updateEsvLambda.

/**
 * Precalculate ESV factors subsequent to lambda propagation.
 */
public void updateEsvLambda() {
    if (!esvTerm) {
        return;
    }
    // Query ExtendedSystem to create local preloads of all lambda quantities.
    numESVs = esvSystem.size();
    if (esvLambda == null || esvLambda.length < nAtoms) {
        esvAtomsScaled = new boolean[nAtoms];
        esvAtomsScaledAlpha = new boolean[nAtoms];
        esvLambda = new double[nAtoms];
        esvIndex = new Integer[nAtoms];
        fill(esvAtomsScaled, false);
        fill(esvAtomsScaledAlpha, false);
        fill(esvLambda, 1.0);
        fill(esvIndex, null);
    }
    /* Preload components for permanent electrostatics. */
    for (int i = 0; i < nAtoms; i++) {
        esvAtomsScaled[i] = esvSystem.isExtended(i);
        esvAtomsScaledAlpha[i] = esvSystem.isAlphaScaled(i);
        esvLambda[i] = esvSystem.getLambda(i);
        esvIndex[i] = (esvSystem.getEsvIndex(i) != null) ? esvSystem.getEsvIndex(i) : null;
    }
    // For atoms with both a foreground and background multipole, preload interpolated multipole and derivative.
    // Initialize unscaled variables.
    esvMultipoleDot = new double[nSymm][nAtoms][10];
    cartMultipoleDotPhi = new double[nAtoms][tensorCount];
    unscaledPolarizability = new double[nAtoms];
    for (int i = 0; i < nAtoms; i++) {
        final Atom ai = atoms[i];
        if (ai.getPolarizeType() == null) {
            logger.warning("Null polarize type during ESV init.");
            continue;
        }
        unscaledPolarizability[i] = ai.getUnscaledPolarizability();
        polarizability[i] = ai.getScaledPolarizability();
    }
}
Also used : Atom(ffx.potential.bonded.Atom)

Example 67 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class ParticleMeshEwaldQI method initSoftCore.

/**
 * Initialize a boolean array of soft atoms and, if requested, ligand vapor
 * electrostatics.
 */
private void initSoftCore(boolean print) {
    if (initSoftCore) {
        return;
    }
    /**
     * Initialize a boolean array that marks soft atoms.
     */
    StringBuilder sb = new StringBuilder("\n Softcore Atoms:\n");
    int count = 0;
    for (int i = 0; i < nAtoms; i++) {
        Atom ai = atoms[i];
        if (ai.applyLambda()) {
            isSoft[i] = true;
            if (print) {
                sb.append(ai.toString()).append("\n");
            }
            count++;
        }
    }
    if (count > 0 && print) {
        logger.info(sb.toString());
    }
    /**
     * Initialize boundary conditions, an n^2 neighbor list and parallel
     * scheduling for ligand vapor electrostatics.
     */
    if (doLigandVaporElec) {
        double maxr = 10.0;
        for (int i = 0; i < nAtoms; i++) {
            Atom ai = atoms[i];
            if (ai.applyLambda()) {
                /**
                 * Determine ligand size.
                 */
                for (int j = i + 1; j < nAtoms; j++) {
                    Atom aj = atoms[j];
                    if (aj.applyLambda()) {
                        double dx = ai.getX() - aj.getX();
                        double dy = ai.getY() - aj.getY();
                        double dz = ai.getZ() - aj.getZ();
                        double r = sqrt(dx * dx + dy * dy + dz * dz);
                        maxr = max(r, maxr);
                    }
                }
            }
        }
        double vacuumOff = 2 * maxr;
        vaporCrystal = new Crystal(3 * vacuumOff, 3 * vacuumOff, 3 * vacuumOff, 90.0, 90.0, 90.0, "P1");
        vaporCrystal.setAperiodic(true);
        NeighborList vacuumNeighborList = new NeighborList(null, vaporCrystal, atoms, vacuumOff, 2.0, parallelTeam);
        vacuumNeighborList.setIntermolecular(false, molecule);
        vaporLists = new int[1][nAtoms][];
        double[][] coords = new double[1][nAtoms * 3];
        for (int i = 0; i < nAtoms; i++) {
            coords[0][i * 3] = atoms[i].getX();
            coords[0][i * 3 + 1] = atoms[i].getY();
            coords[0][i * 3 + 2] = atoms[i].getZ();
        }
        vacuumNeighborList.buildList(coords, vaporLists, isSoft, true, true);
        vaporPermanentSchedule = vacuumNeighborList.getPairwiseSchedule();
        vaporEwaldSchedule = vaporPermanentSchedule;
        vacuumRanges = new Range[maxThreads];
        vacuumNeighborList.setDisableUpdates(forceField.getBoolean(ForceField.ForceFieldBoolean.DISABLE_NEIGHBOR_UPDATES, false));
    } else {
        vaporCrystal = null;
        vaporLists = null;
        vaporPermanentSchedule = null;
        vaporEwaldSchedule = null;
        vacuumRanges = null;
    }
    /**
     * Set this flag to true to avoid re-initialization.
     */
    initSoftCore = true;
}
Also used : Atom(ffx.potential.bonded.Atom) Crystal(ffx.crystal.Crystal)

Example 68 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class ParticleMeshEwaldQI method initAtomArrays.

private void initAtomArrays() {
    if (localMultipole == null || localMultipole.length < nAtoms || lambdaTerm || esvTerm) {
        localMultipole = new double[nAtoms][10];
        frame = new MultipoleType.MultipoleFrameDefinition[nAtoms];
        axisAtom = new int[nAtoms][];
        cartMultipolePhi = new double[nAtoms][tensorCount];
        directDipole = new double[nAtoms][3];
        directDipoleCR = new double[nAtoms][3];
        cartesianDipolePhi = new double[nAtoms][tensorCount];
        cartesianDipolePhiCR = new double[nAtoms][tensorCount];
        ip11 = new int[nAtoms][];
        ip12 = new int[nAtoms][];
        ip13 = new int[nAtoms][];
        thole = new double[nAtoms];
        ipdamp = new double[nAtoms];
        polarizability = new double[nAtoms];
        realSpaceSchedule = new PairwiseSchedule(maxThreads, nAtoms, realSpaceRanges);
        if (scfAlgorithm == SCFAlgorithm.CG) {
            rsd = new double[3][nAtoms];
            rsdCR = new double[3][nAtoms];
            rsdPre = new double[3][nAtoms];
            rsdPreCR = new double[3][nAtoms];
            conj = new double[3][nAtoms];
            conjCR = new double[3][nAtoms];
            vec = new double[3][nAtoms];
            vecCR = new double[3][nAtoms];
        }
        /**
         * Initialize per-thread memory for collecting the gradient, torque,
         * field and chain-rule field.
         */
        grad = new double[maxThreads][3][nAtoms];
        torque = new double[maxThreads][3][nAtoms];
        field = new double[maxThreads][3][nAtoms];
        fieldCR = new double[maxThreads][3][nAtoms];
        if (lambdaTerm) {
            lambdaGrad = new double[maxThreads][3][nAtoms];
            lambdaTorque = new double[maxThreads][3][nAtoms];
        }
        isSoft = new boolean[nAtoms];
        fill(isSoft, false);
        use = new boolean[nAtoms];
        fill(use, true);
        coordinates = new double[nSymm][3][nAtoms];
        globalMultipole = new double[nSymm][nAtoms][10];
        inducedDipole = new double[nSymm][nAtoms][3];
        inducedDipoleCR = new double[nSymm][nAtoms][3];
        if (scfPredictor != null) {
            if (esvTerm) {
                throw new UnsupportedOperationException();
            }
            scfPredictor.setInducedDipoleReferences(inducedDipole, inducedDipoleCR, lambdaTerm);
        }
        /* ESV flag array, initialized regardless of esvTerm. */
        // True for other ESV residue atoms.
        esvAtomsScaled = new boolean[nAtoms];
        esvAtomsScaledAlpha = new boolean[nAtoms];
        fill(esvAtomsScaled, false);
        fill(esvAtomsScaledAlpha, false);
        lambdaFactors = new LambdaFactors[maxThreads];
        for (int i = 0; i < maxThreads; i++) {
            if (esvTerm) {
                // Invoked every time through inner loops.
                lambdaFactors[i] = new LambdaFactorsESV();
            } else if (lambdaTerm) {
                // Invoked on calls to setLambda().
                lambdaFactors[i] = new LambdaFactorsOSRW();
            } else {
                // Invoked never; inoperative defaults.
                lambdaFactors[i] = LambdaDefaults;
            }
        }
        /**
         * The size of reduced neighbor list depends on the size of the real
         * space cutoff.
         */
        realSpaceLists = new int[nSymm][nAtoms][];
        realSpaceCounts = new int[nSymm][nAtoms];
        preconditionerLists = new int[nSymm][nAtoms][preconditionerListSize];
        preconditionerCounts = new int[nSymm][nAtoms];
    }
    /**
     * Assign multipole parameters and polarization groups.
     */
    for (int i = 0; i < nAtoms; i++) {
        MultipoleType.multipoleTypeFactory(atoms[i], forceField);
        localMultipole[i] = atoms[i].getMultipoleType().getMultipole();
        axisAtom[i] = atoms[i].getAxisAtomIndices();
        frame[i] = atoms[i].getMultipoleType().frameDefinition;
    }
    /**
     * Assign polarization groups.
     */
    assignPolarizationGroups();
    /**
     * Fill the thole, inverse polarization damping and polarizability
     * arrays.
     */
    for (Atom ai : atoms) {
        PolarizeType polarizeType = ai.getPolarizeType();
        int index = ai.getIndex() - 1;
        thole[index] = polarizeType.thole;
        ipdamp[index] = polarizeType.pdamp;
        if (!(ipdamp[index] > 0.0)) {
            ipdamp[index] = Double.POSITIVE_INFINITY;
        } else {
            ipdamp[index] = 1.0 / ipdamp[index];
        }
        polarizability[index] = polarizeType.polarizability;
    }
    if (esvTerm) {
        updateEsvLambda();
    }
}
Also used : PolarizeType(ffx.potential.parameters.PolarizeType) MultipoleType(ffx.potential.parameters.MultipoleType) Atom(ffx.potential.bonded.Atom)

Example 69 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class ExtendedSystem method addVariable.

/**
 * Prefer ExtendedSystem::populate to manual ESV creation.
 */
public void addVariable(ExtendedVariable esv) {
    esvSystemActive = true;
    if (esvList == null) {
        esvList = new ArrayList<>();
    }
    if (esvList.contains(esv)) {
        logger.warning(format("Attempted to add duplicate variable %s to system.", esv.toString()));
        return;
    }
    esvList.add(esv);
    numESVs = esvList.size();
    if (esv instanceof TitrationESV) {
        if (constantSystemPh == null) {
            logger.severe("Set ExtendedSystem (constant) pH before adding TitrationESVs.");
        }
        phTerm = true;
    }
    for (int i = 0; i < nAtomsExt; i++) {
        if (esv.viewUnsharedAtoms().contains(extendedAtoms[i])) {
            esvForUnshared[i] = esv;
        } else if (esv.viewSharedAtoms().contains(extendedAtoms[i])) {
            esvForShared[i] = esv;
        }
    }
    fg2bgIdx = new int[nAtomsExt];
    for (int i = 0; i < nAtomsExt; i++) {
        Atom atom = extendedAtoms[i];
        if (isExtended(i)) {
            Atom bg = atom.getEsv().getBackgroundForAtom(atom);
            fg2bgIdx[i] = (bg != null) ? bg.getIndex() : -1;
        }
    }
    updateListeners();
}
Also used : Atom(ffx.potential.bonded.Atom)

Example 70 with Atom

use of ffx.potential.bonded.Atom in project ffx by mjschnie.

the class TitrationUtils method propagateInactiveResidues.

/**
 * Copies atomic coordinates from each active residue to its inactive
 * counterparts. Inactive hydrogen coordinates are updated by geometry with
 * the propagated heavies.
 */
public static void propagateInactiveResidues(MultiResidue multiRes, boolean propagateDynamics) {
    // Propagate all atom coordinates from active residues to their inactive counterparts.
    Residue active = multiRes.getActive();
    String activeResName = active.getName();
    List<Residue> inactives = multiRes.getInactive();
    for (Atom activeAtom : active.getAtomList()) {
        String activeName = activeAtom.getName();
        for (Residue inactive : inactives) {
            Atom inactiveAtom = (Atom) inactive.getAtomNode(activeName);
            if (inactiveAtom != null) {
                // Propagate position and gradient.
                double[] activeXYZ = activeAtom.getXYZ(null);
                inactiveAtom.setXYZ(activeXYZ);
                double[] grad = new double[3];
                activeAtom.getXYZGradient(grad);
                inactiveAtom.setXYZGradient(grad[0], grad[1], grad[2]);
                if (propagateDynamics) {
                    // Propagate velocity, acceleration, and previous acceleration.
                    double[] activeVelocity = new double[3];
                    activeAtom.getVelocity(activeVelocity);
                    inactiveAtom.setVelocity(activeVelocity);
                    double[] activeAccel = new double[3];
                    activeAtom.getAcceleration(activeAccel);
                    inactiveAtom.setAcceleration(activeAccel);
                    double[] activePrevAcc = new double[3];
                    activeAtom.getPreviousAcceleration(activePrevAcc);
                    inactiveAtom.setPreviousAcceleration(activePrevAcc);
                }
            } else {
                if (activeName.equals("C") || activeName.equals("O") || activeName.equals("N") || activeName.equals("CA") || activeName.equals("H") || activeName.equals("HA")) {
                // Backbone atoms aren't supposed to exist in inactive multiResidue components; so no problem.
                } else if (isTitratableHydrogen(activeAtom)) {
                /**
                 * i.e. ((activeResName.equals("LYS") &&
                 * activeName.equals("HZ3")) ||
                 * (activeResName.equals("TYR") &&
                 * activeName.equals("HH")) ||
                 * (activeResName.equals("CYS") &&
                 * activeName.equals("HG")) ||
                 * (activeResName.equals("HIS") &&
                 * (activeName.equals("HD1") ||
                 * activeName.equals("HE2"))) ||
                 * (activeResName.equals("HID") &&
                 * activeName.equals("HD1")) ||
                 * (activeResName.equals("HIE") &&
                 * activeName.equals("HE2")) ||
                 * (activeResName.equals("ASH") &&
                 * activeName.equals("HD2")) ||
                 * (activeResName.equals("GLH") &&
                 * activeName.equals("HE2")))
                 */
                // These titratable protons are handled below; so no problem.
                } else {
                    // Now we have a problem.
                    logger.warning(format("Couldn't propagate inactive MultiResidue atom: %s: %s, %s", multiRes, activeName, activeAtom));
                }
            }
        }
    }
    rebuildStrandedProtons(multiRes);
}
Also used : MultiResidue(ffx.potential.bonded.MultiResidue) Residue(ffx.potential.bonded.Residue) Atom(ffx.potential.bonded.Atom)

Aggregations

Atom (ffx.potential.bonded.Atom)206 Residue (ffx.potential.bonded.Residue)42 Bond (ffx.potential.bonded.Bond)37 CoordRestraint (ffx.potential.nonbonded.CoordRestraint)34 ArrayList (java.util.ArrayList)33 OpenMM_System_addConstraint (simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint)23 Polymer (ffx.potential.bonded.Polymer)22 IOException (java.io.IOException)22 MSNode (ffx.potential.bonded.MSNode)21 Crystal (ffx.crystal.Crystal)20 Molecule (ffx.potential.bonded.Molecule)19 File (java.io.File)19 MultiResidue (ffx.potential.bonded.MultiResidue)17 MultipoleType (ffx.potential.parameters.MultipoleType)13 MissingHeavyAtomException (ffx.potential.bonded.BondedUtils.MissingHeavyAtomException)12 PointerByReference (com.sun.jna.ptr.PointerByReference)11 MissingAtomTypeException (ffx.potential.bonded.BondedUtils.MissingAtomTypeException)11 AtomType (ffx.potential.parameters.AtomType)11 List (java.util.List)11 MolecularAssembly (ffx.potential.MolecularAssembly)10