<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1">
<style type="text/css" style="display:none;"> P {margin-top:0;margin-bottom:0;} </style>
</head>
<body dir="ltr">
<div style="font-family: Calibri, Arial, Helvetica, sans-serif; font-size: 12pt; color: rgb(0, 0, 0);">
Hi&nbsp;<span style="color:rgb(32, 31, 30);font-family:&quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif;font-size:15px;background-color:rgb(255, 255, 255);display:inline !important">Lorién</span></div>
<div style="font-family: Calibri, Arial, Helvetica, sans-serif; font-size: 12pt; color: rgb(0, 0, 0);">
<span style="color:rgb(32, 31, 30);font-family:&quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif;font-size:15px;background-color:rgb(255, 255, 255);display:inline !important"><br>
</span></div>
<div style="font-family: Calibri, Arial, Helvetica, sans-serif; font-size: 12pt; color: rgb(0, 0, 0);">
<span style="color:rgb(32, 31, 30);font-family:&quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif;font-size:15px;background-color:rgb(255, 255, 255);display:inline !important">One
 the possible solutions (there are atomic numbers stored in the gmx_mtop_t):</span></div>
<div style="font-family: Calibri, Arial, Helvetica, sans-serif; font-size: 12pt; color: rgb(0, 0, 0);">
<br>
</div>
<div style="font-family: Calibri, Arial, Helvetica, sans-serif; font-size: 12pt; color: rgb(0, 0, 0);">
<br>
</div>
<div style=""><font color="#201f1e"><span style="font-size: 15px;">&nbsp; &nbsp; std::vector&lt;int&gt;&nbsp;<span style="background-color:rgb(255, 255, 255);display:inline !important">atomNumbers</span></span></font></div>
<div style=""><span style="color: rgb(32, 31, 30); font-family: &quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif; font-size: 15px; background-color: rgb(255, 255, 255); display: inline !important;">
<div><br>
</div>
<div>&nbsp; &nbsp; AtomIterator atoms(*mtop);</div>
<div>&nbsp; &nbsp; while (atoms-&gt;globalAtomNumber() &lt; mtop-&gt;natoms)</div>
<div>&nbsp; &nbsp; {</div>
<div>&nbsp; &nbsp; &nbsp; &nbsp; atomNumbers.push_back(atoms-&gt;atom().atomnumber);</div>
<div>&nbsp; &nbsp; &nbsp; &nbsp; atoms++;</div>
&nbsp; &nbsp; }<br>
</span></div>
<div style=""><span style="color: rgb(32, 31, 30); font-family: &quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif; font-size: 15px; background-color: rgb(255, 255, 255); display: inline !important;"><br>
</span></div>
<div style=""><span style="color: rgb(32, 31, 30); font-family: &quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif; font-size: 15px; background-color: rgb(255, 255, 255); display: inline !important;">for
 that you need some functionality from&nbsp;&quot;gromacs/topology/mtop_util.h&quot;</span></div>
<div style=""><span style="color: rgb(32, 31, 30); font-family: &quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif; font-size: 15px;">And mtop here is a pointer to gmx_mtop_t structure</span><br>
</div>
<div style=""><br>
</div>
<div style=""><span style="color: rgb(32, 31, 30); font-family: &quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif; font-size: 15px; background-color: rgb(255, 255, 255); display: inline !important;"><br>
</span></div>
<div style=""><span style="color: rgb(32, 31, 30); font-family: &quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif; font-size: 15px; background-color: rgb(255, 255, 255); display: inline !important;">Best,</span></div>
<div style=""><span style="color: rgb(32, 31, 30); font-family: &quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif; font-size: 15px; background-color: rgb(255, 255, 255); display: inline !important;">Dmitry</span></div>
<div style="font-family: Calibri, Arial, Helvetica, sans-serif; font-size: 12pt; color: rgb(0, 0, 0);">
<span style="color:rgb(32, 31, 30);font-family:&quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif;font-size:15px;background-color:rgb(255, 255, 255);display:inline !important"><br>
</span></div>
<div style="font-family: Calibri, Arial, Helvetica, sans-serif; font-size: 12pt; color: rgb(0, 0, 0);">
<span style="color:rgb(32, 31, 30);font-family:&quot;Segoe UI&quot;, &quot;Segoe UI Web (West European)&quot;, &quot;Segoe UI&quot;, -apple-system, BlinkMacSystemFont, Roboto, &quot;Helvetica Neue&quot;, sans-serif;font-size:15px;background-color:rgb(255, 255, 255);display:inline !important"><br>
</span></div>
<div id="appendonsend"></div>
<hr style="display:inline-block;width:98%" tabindex="-1">
<div id="divRplyFwdMsg" dir="ltr"><font face="Calibri, sans-serif" style="font-size:11pt" color="#000000"><b>From:</b> gromacs.org_gmx-developers-bounces@maillist.sys.kth.se &lt;gromacs.org_gmx-developers-bounces@maillist.sys.kth.se&gt; on behalf of Lorién López
 Villellas &lt;lorien.lopez@bsc.es&gt;<br>
<b>Sent:</b> Thursday, February 11, 2021 17:22<br>
<b>To:</b> gmx-developers@gromacs.org &lt;gmx-developers@gromacs.org&gt;<br>
<b>Cc:</b> Pablo Ibañez &lt;imarin@unizar.es&gt;; Miquel Moreto &lt;miquel.moreto@bsc.es&gt;; SANTIAGO MARCO SOLA &lt;santiago.marco@bsc.es&gt;; Jesús Alastruey &lt;jalastru@unizar.es&gt;; Kanduri Prashanth &lt;pkanduri@cscs.ch&gt;; gromacs.org_gmx-developers@maillist.sys.kth.se &lt;gromacs.org_gmx-developers@maillist.sys.kth.se&gt;;
 Keller Sebastian &lt;keller@cscs.ch&gt;; Holanda Rusu Victor &lt;victor.holanda@cscs.ch&gt;<br>
<b>Subject:</b> Re: [gmx-developers] How to obtain the element of an atom?</font>
<div>&nbsp;</div>
</div>
<div>
<div dir="ltr">Hi.<br>
<br>
<blockquote class="x_gmail_quote" style="margin:0px 0px 0px 0.8ex; border-left:1px solid rgb(204,204,204); padding-left:1ex">
Why is the element number relevant for constraints? The element number is never used in MD calculations.</blockquote>
<br>
We have performed a heavy code optimization based on the structural patterns of the peptide chains. In order to apply this optimization, we need a very specific bond numbering. To get to this bond numbering, we use a bond reordering phase. The reordering algorithm
 needs to know the elements of the molecule's atoms to locate itself.<br>
<br>
<blockquote class="x_gmail_quote" style="margin:0px 0px 0px 0.8ex; border-left:1px solid rgb(204,204,204); padding-left:1ex">
Working with the legacy gromacs topology data structures, such as the iatoms (really just a vector of ints) and iparams can be quite challenging. Additionally, adding new code to constraints code is likely to face challenges based on the core gromacs team being
 hesitant to take on new maintenance burdens. However, depending on how much effort you have to spend on this, I might suggest that you could check out the nblib topology and force calculation API in gromacs 2021. Using this, it could be much easier to implement
 your new algorithm. Please have a look at the sample scripts in api/nblib/samples to get an idea of how, for instance, listed forces such as bonds are handled by nblib. If this looks promising, perhaps we can have a chat about potential integration pathways.</blockquote>
<br>
This project has been running for some years now (since 2013), so we started working with the old C-GROMACS data structures. So far, we have used a custom benchmark, but now we want to evaluate our algorithm's performance in a real environment. I will take
 a look at nblib.<br>
<br>
Regards,<br>
Lorién<br>
</div>
<br>
<div class="x_gmail_quote">
<div dir="ltr" class="x_gmail_attr">El jue, 11 feb 2021 a las 15:52, Joe Jordan (&lt;<a href="mailto:e.jjordan12@gmail.com">e.jjordan12@gmail.com</a>&gt;) escribió:<br>
</div>
<blockquote class="x_gmail_quote" style="margin:0px 0px 0px 0.8ex; border-left:1px solid rgb(204,204,204); padding-left:1ex">
<div dir="ltr">Hi Lorién,
<div><br>
</div>
<div>Working with the legacy gromacs topology data structures, such as the iatoms (really just a vector of ints) and iparams can be quite challenging. Additionally, adding new code to constraints code is likely to face challenges based on the core gromacs team
 being hesitant to take on new maintenance burdens. However, depending on how much effort you have to spend on this, I might suggest that you could check out the nblib topology and force calculation API in gromacs 2021. Using this, it could be much easier to
 implement your new algorithm. Please have a look at the sample scripts in api/nblib/samples to get an idea of how, for instance, listed forces such as bonds are handled by nblib. If this looks promising, perhaps we can have a chat about potential integration
 pathways.</div>
<div><br>
</div>
<div>Joe</div>
</div>
<br>
<div class="x_gmail_quote">
<div dir="ltr" class="x_gmail_attr">On Thu, Feb 11, 2021 at 3:30 PM Lorién López Villellas &lt;<a href="mailto:lorien.lopez@bsc.es" target="_blank">lorien.lopez@bsc.es</a>&gt; wrote:<br>
</div>
<blockquote class="x_gmail_quote" style="margin:0px 0px 0px 0.8ex; border-left:1px solid rgb(204,204,204); padding-left:1ex">
<div dir="ltr">Hi all.<br>
<br>
We are working on a new parallel algorithm for imposing constraints based on SHAKE. We are currently trying to integrate it into GROMACS.<br>
<br>
We use the iatoms array (topology of triplets of constraint type, and indices of the two atoms that compose a bond) just like SHAKE. We need a way to obtain the element (Carbon, Nitrogen, Hidrogen...) of each of the atom indexes in iatoms. As far as we have
 researched, there is no way to access that data from src/gromacs/mdlib/constr. Where is that information stored?<br>
<br>
If there is no way to obtain an atom element, we could also adapt our code to work with similar information if it exists.<br>
<br>
We really appreciate any help you can provide.&nbsp;<br>
</div>
-- <br>
Gromacs Developers mailing list<br>
<br>
* Please search the archive at <a href="http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_List" rel="noreferrer" target="_blank">
http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_List</a> before posting!<br>
<br>
* Can't post? Read <a href="http://www.gromacs.org/Support/Mailing_Lists" rel="noreferrer" target="_blank">
http://www.gromacs.org/Support/Mailing_Lists</a><br>
<br>
* For (un)subscribe requests visit<br>
<a href="https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-developers" rel="noreferrer" target="_blank">https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-developers</a> or send a mail to
<a href="mailto:gmx-developers-request@gromacs.org" target="_blank">gmx-developers-request@gromacs.org</a>.</blockquote>
</div>
-- <br>
Gromacs Developers mailing list<br>
<br>
* Please search the archive at <a href="http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_List" rel="noreferrer" target="_blank">
http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_List</a> before posting!<br>
<br>
* Can't post? Read <a href="http://www.gromacs.org/Support/Mailing_Lists" rel="noreferrer" target="_blank">
http://www.gromacs.org/Support/Mailing_Lists</a><br>
<br>
* For (un)subscribe requests visit<br>
<a href="https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-developers" rel="noreferrer" target="_blank">https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-developers</a> or send a mail to
<a href="mailto:gmx-developers-request@gromacs.org" target="_blank">gmx-developers-request@gromacs.org</a>.</blockquote>
</div>
</div>
</body>
</html>