<div dir="ltr"><div>Jan,</div><div><br></div><div>If you want to profile and try to tune performance, I suggest you focus on the main simulation tool called "mdrun" and in particular its core functionalities as this is the tool that will be consuming the vast majority of compute cycles.<br></div><div><br></div><div>If you are using profiling tools to explore the code, find hot functions and such, I suggest that you grab one of the benchmark systems from here:</div><div><a href="http://ftp.gromacs.org/benchmarks/">http://ftp.gromacs.org/benchmarks/</a></div><div>and run one of those (make sure to pick a "PME" input if there are multiple).<br></div><div><br></div><div>And start by running those through a profiler ( e.g. run "gmx mdrun -s topol.tpr -nsteps 1000" for a 1000-iteration run). If you are proficient at profiling tools, you can also drop in manual annotation (or profiler start/stop) into out internal timer/profiler facilities (see src/gromacs/timing/wallcycle.cpp).</div><div><br></div><div>Let us know if you have further questions.<br></div><div><br></div><div>Cheers,<br></div><div><div><div dir="ltr" class="gmail_signature" data-smartmail="gmail_signature">--<br>Szilárd</div></div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, Mar 23, 2020 at 3:36 PM jan <<a href="mailto:rtm443x@googlemail.com">rtm443x@googlemail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">Hi, see below<br>
<br>
On 23/03/2020, cblau <<a href="mailto:cblau@gwdg.de" target="_blank">cblau@gwdg.de</a>> wrote:<br>
<br>
[snip]<br>
><br>
> As GROMACS is already highly specialised, it is a bit challenging to<br>
> speed up quickly (though definitely possible to contribute in meaningful<br>
> ways).<br>
<br>
My initial look is profile then see look at the most basic things such<br>
as memory hierarchy use. I mean, really basic stuff which you'll<br>
already have covered but at least it's a learning phase for me.<br>
<br>
><br>
> With all urgency in mind, there is quite a time-delay in the way new<br>
> code gets incorporated, because we do very much require scientific<br>
> correctness and thorough code review to make sure that all simulation<br>
> results are meaningful and we don't waste a ton of people efforts and<br>
> compute power on a slight bug.<br>
<br>
Is good<br>
<br>
><br>
> The largest speed up in a short time frame (the next months) is inmy<br>
> option achieved by teaching people how to use GROMACS in the most<br>
> efficient ways and make sure that people running simulations run them<br>
> correctly.<br>
<br>
Knowing how to use a tool correctly is best way to get good results. I<br>
keep seeing the results where this does not happen. So I get you.<br>
I'll look on the site for reading material but any other pointers<br>
would help - but bear in mind I'm very new to this, and numerical<br>
simulation isn't my thing so start simple!<br>
<br>
><br>
> If you like, we can briefly chat tomorrow on how GROMACS is set up and<br>
> the "usual ways" to contribute to GROMACS if you send me a personal email.<br>
><br>
<br>
Will do<br>
<br>
thanks<br>
<br>
jan<br>
<br>
<br>
><br>
> Best,<br>
> Christian<br>
><br>
> On 2020-03-23 14:58, jan wrote:<br>
>> Hi,<br>
>> I'm a general back-end dev. Given the situation, and folding@home<br>
>> using gromacs, I thought I'd poke through the code. I noticed<br>
>> something unexpected, and was advised to email it here. in edsam.cpp,<br>
>> this:<br>
>><br>
>><br>
>> void do_linacc(rvec* xcoll, t_edpar* edi)<br>
>> {<br>
>> /* loop over linacc vectors */<br>
>> for (int i = 0; i < edi->vecs.linacc.neig; i++)<br>
>> {<br>
>> /* calculate the projection */<br>
>> real proj = projectx(*edi, xcoll, edi->vecs.linacc.vec[i]);<br>
>><br>
>><br>
>> /* calculate the correction */<br>
>> real preFactor = 0.0;<br>
>> if (edi->vecs.linacc.stpsz[i] > 0.0)<br>
>> {<br>
>> if ((proj - edi->vecs.linacc.refproj[i]) < 0.0)<br>
>> {<br>
>> preFactor = edi->vecs.linacc.refproj[i] - proj;<br>
>> }<br>
>> }<br>
>> if (edi->vecs.linacc.stpsz[i] < 0.0)<br>
>> {<br>
>> if ((proj - edi->vecs.linacc.refproj[i]) > 0.0)<br>
>> {<br>
>> preFactor = edi->vecs.linacc.refproj[i] - proj;<br>
>> }<br>
>> }<br>
>> [...]<br>
>><br>
>><br>
>> In both cases it reaches the same code<br>
>><br>
>> preFactor = edi->vecs.linacc.refproj[i] - proj<br>
>><br>
>> That surprised me a bit, is it deliberate? If so it may be the code<br>
>> can be simplified anyway.<br>
>><br>
>> That aside, if you're looking for performance I might be able to help.<br>
>> I don't know the high level stuff *at this point* and my C++ is so<br>
>> rusty it creaks, but I can brush that up, do profiling and whatnot.<br>
>> I'm pretty experience, just not in this area. Speeding things up is<br>
>> something I've got a track record of (though I usually have a good<br>
>> feel for the problem domain first, which I don't here)<br>
>><br>
>> Would it be of some value for me to try getting more speed? If so,<br>
>> first thing I'd need is to get this running under cygwin, which I'm<br>
>> struggling with.<br>
>><br>
>> regards<br>
>><br>
>> jan<br>
>><br>
> --<br>
> Gromacs Developers mailing list<br>
><br>
> * Please search the archive at<br>
> <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<br>
> 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<br>
> send a mail to <a href="mailto:gmx-developers-request@gromacs.org" target="_blank">gmx-developers-request@gromacs.org</a>.<br>
><br>
-- <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>.<br>
</blockquote></div>