Minima_4.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "16b8ce47-81e0-49c9-9e52-fd3d6a12fc73",
   "metadata": {},
   "outputs": [],
   "source": [
    "from openmm.app import *\n",
    "from openmm import *\n",
    "from openmm import unit as u\n",
    "import sys\n",
    "import mdtraj as md\n",
    "import numpy as np\n",
    "from tqdm import tqdm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "2cba4651-444f-44bd-a569-9ab4aa9ab555",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Nombre del archivo PDB del péptido\n",
    "pdb_filename = \"pdbs/V5.pdb\"\n",
    "run_index = 0\n",
    "steps_segment = 1000\n",
    "n_segments = 1000\n",
    "temperature = 500 * u.kelvin"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "53596979-342c-459c-82f5-d384987d4ee4",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Cargar la estructura del péptido desde el archivo PDB\n",
    "pdb = PDBFile(pdb_filename)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "54e5f0a3-de71-4813-83d4-141b70aea31a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Definir la fuerza de campo (AMBER ff14SB es común para péptidos)\n",
    "forcefield = ForceField(\"amber14-all.xml\", \"amber14/tip3p.xml\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "3cebe2fc-c3d6-4947-b80a-3720c1d84001",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Crear el sistema en vacío (sin agua ni iones)\n",
    "system = forcefield.createSystem(\n",
    "    pdb.topology,\n",
    "    nonbondedMethod=NoCutoff,  # Sin tratamiento de fronteras periódicas\n",
    "    constraints=HBonds         # Restringe enlaces de hidrógeno para un paso de integración más largo\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "82d23f65-9107-4ceb-b55c-4aa9b8e0ff0f",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Configurar el integrador (Langevin dynamics)\n",
    "friction = 1 / u.picosecond\n",
    "timestep = 2 * u.femtoseconds\n",
    "integrator = LangevinIntegrator(temperature, friction, timestep)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "ef2c24dd-c3d6-49d2-b552-62f60b06cb1e",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Crear una plataforma para ejecutar la simulación (e.g., CUDA si tienes GPU)\n",
    "platform = Platform.getPlatformByName(\"CUDA\")  # Usa \"CPU\" si no tienes GPU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "4661c388-e453-48f5-b164-db907b5b635c",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Configurar la simulación\n",
    "simulation = Simulation(pdb.topology, system, integrator, platform)\n",
    "simulation.context.setPositions(pdb.positions)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "8cc726bb-234c-4b9c-bb63-fc296d1ade6e",
   "metadata": {},
   "outputs": [],
   "source": [
    "state = simulation.context.getState(getEnergy=True)\n",
    "energy = state.getPotentialEnergy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "11d5ced6-1f6b-49b9-b8f9-400cd9eda7df",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Minimización de energía\n",
    "simulation.minimizeEnergy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "9355ce22-0f03-4b89-a1ec-466580945c42",
   "metadata": {},
   "outputs": [],
   "source": [
    "md_topology = md.Topology.from_openmm(simulation.topology)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "3cffdef3-fdfb-4146-91e8-88a5296f1785",
   "metadata": {},
   "outputs": [],
   "source": [
    "def es_el_mismo(angs1, angs2):\n",
    "    if np.max(np.abs(angs1-angs2))<0.1:\n",
    "        return True\n",
    "    else:\n",
    "        return False"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "e20d2340-ab9d-4cdd-b29a-9352b987405a",
   "metadata": {},
   "outputs": [],
   "source": [
    "traj_inh = []\n",
    "dihed_db = []\n",
    "coors_db = []\n",
    "energy_db = []"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "064235bf-a800-4f84-88a3-194d10ff2d91",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 1000/1000 [00:58<00:00, 17.15it/s]\n"
     ]
    }
   ],
   "source": [
    "for ii in tqdm(range(n_segments)):\n",
    "    simulation.step(steps_segment)\n",
    "    state = simulation.context.getState(getEnergy=True)\n",
    "    energy = state.getPotentialEnergy()\n",
    "    simulation.minimizeEnergy()\n",
    "    min_state = simulation.context.getState(getEnergy=True, getPositions=True)\n",
    "    min_energy = min_state.getPotentialEnergy()\n",
    "    min_positions = min_state.getPositions(asNumpy=True)\n",
    "    traj = md.Trajectory(min_positions / u.nanometer, md_topology)\n",
    "    phis = md.compute_phi(traj)[1]\n",
    "    psis = md.compute_psi(traj)[1]\n",
    "    dihed_angs=np.concatenate((phis[0],psis[0]))\n",
    "    visitado = False\n",
    "    for unique_index in range(len(dihed_db)):\n",
    "        aux = es_el_mismo(dihed_angs, dihed_db[unique_index])\n",
    "        if aux == True:\n",
    "            visitado = True\n",
    "            traj_inh.append(unique_index)\n",
    "            break\n",
    "    if visitado == False:\n",
    "        traj_inh.append(len(dihed_db))\n",
    "        dihed_db.append(dihed_angs)\n",
    "        coors_db.append(min_positions._value)\n",
    "        energy_db.append(min_energy._value)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "f76fec26-2a7b-476b-85ff-4c85c10b320d",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pickle as pickle\n",
    "import gzip"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "195e741a-650a-46d3-9eb4-4af9c2458a8d",
   "metadata": {},
   "outputs": [],
   "source": [
    "with gzip.open(f'energies_{run_index}.pkl.gz', 'wb') as fff:\n",
    "    pickle.dump(energy_db, fff)\n",
    "\n",
    "with gzip.open(f'traj_inh_{run_index}.pkl.gz', 'wb') as fff:\n",
    "    pickle.dump(traj_inh, fff)\n",
    "\n",
    "with gzip.open(f'coors_{run_index}.pkl.gz', 'wb') as fff:\n",
    "    pickle.dump(coors_db, fff)\n",
    "\n",
    "with gzip.open(f'dihed_{run_index}.pkl.gz', 'wb') as fff:\n",
    "    pickle.dump(dihed_db, fff)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "ccec410a-ee3f-48af-b7ed-8b5a8b3c4b7c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "960000"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "60000*16"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "1303f84f-4282-412e-8197-e3f1fb8335e0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2000000"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "1000000*2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "aee7ee4d-dc09-430a-9eee-0cb8e2d00589",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}