1
0
mirror of https://github.com/TREX-CoE/qmckl.git synced 2024-07-22 18:57:40 +02:00
qmckl/qmckl_blas.html

863 lines
35 KiB
HTML
Raw Normal View History

<?xml version="1.0" encoding="utf-8"?>
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Strict//EN"
"http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd">
<html xmlns="http://www.w3.org/1999/xhtml" lang="en" xml:lang="en">
<head>
<!-- 2022-01-08 Sat 14:20 -->
<meta http-equiv="Content-Type" content="text/html;charset=utf-8" />
<meta name="viewport" content="width=device-width, initial-scale=1" />
<title>BLAS functions</title>
<meta name="generator" content="Org mode" />
<meta name="author" content="TREX CoE" />
<style type="text/css">
<!--/*--><![CDATA[/*><!--*/
.title { text-align: center;
margin-bottom: .2em; }
.subtitle { text-align: center;
font-size: medium;
font-weight: bold;
margin-top:0; }
.todo { font-family: monospace; color: red; }
.done { font-family: monospace; color: green; }
.priority { font-family: monospace; color: orange; }
.tag { background-color: #eee; font-family: monospace;
padding: 2px; font-size: 80%; font-weight: normal; }
.timestamp { color: #bebebe; }
.timestamp-kwd { color: #5f9ea0; }
.org-right { margin-left: auto; margin-right: 0px; text-align: right; }
.org-left { margin-left: 0px; margin-right: auto; text-align: left; }
.org-center { margin-left: auto; margin-right: auto; text-align: center; }
.underline { text-decoration: underline; }
#postamble p, #preamble p { font-size: 90%; margin: .2em; }
p.verse { margin-left: 3%; }
pre {
border: 1px solid #ccc;
box-shadow: 3px 3px 3px #eee;
padding: 8pt;
font-family: monospace;
overflow: auto;
margin: 1.2em;
}
pre.src {
position: relative;
overflow: visible;
padding-top: 1.2em;
}
pre.src:before {
display: none;
position: absolute;
background-color: white;
top: -10px;
right: 10px;
padding: 3px;
border: 1px solid black;
}
pre.src:hover:before { display: inline;}
/* Languages per Org manual */
pre.src-asymptote:before { content: 'Asymptote'; }
pre.src-awk:before { content: 'Awk'; }
pre.src-C:before { content: 'C'; }
/* pre.src-C++ doesn't work in CSS */
pre.src-clojure:before { content: 'Clojure'; }
pre.src-css:before { content: 'CSS'; }
pre.src-D:before { content: 'D'; }
pre.src-ditaa:before { content: 'ditaa'; }
pre.src-dot:before { content: 'Graphviz'; }
pre.src-calc:before { content: 'Emacs Calc'; }
pre.src-emacs-lisp:before { content: 'Emacs Lisp'; }
pre.src-fortran:before { content: 'Fortran'; }
pre.src-gnuplot:before { content: 'gnuplot'; }
pre.src-haskell:before { content: 'Haskell'; }
pre.src-hledger:before { content: 'hledger'; }
pre.src-java:before { content: 'Java'; }
pre.src-js:before { content: 'Javascript'; }
pre.src-latex:before { content: 'LaTeX'; }
pre.src-ledger:before { content: 'Ledger'; }
pre.src-lisp:before { content: 'Lisp'; }
pre.src-lilypond:before { content: 'Lilypond'; }
pre.src-lua:before { content: 'Lua'; }
pre.src-matlab:before { content: 'MATLAB'; }
pre.src-mscgen:before { content: 'Mscgen'; }
pre.src-ocaml:before { content: 'Objective Caml'; }
pre.src-octave:before { content: 'Octave'; }
pre.src-org:before { content: 'Org mode'; }
pre.src-oz:before { content: 'OZ'; }
pre.src-plantuml:before { content: 'Plantuml'; }
pre.src-processing:before { content: 'Processing.js'; }
pre.src-python:before { content: 'Python'; }
pre.src-R:before { content: 'R'; }
pre.src-ruby:before { content: 'Ruby'; }
pre.src-sass:before { content: 'Sass'; }
pre.src-scheme:before { content: 'Scheme'; }
pre.src-screen:before { content: 'Gnu Screen'; }
pre.src-sed:before { content: 'Sed'; }
pre.src-sh:before { content: 'shell'; }
pre.src-sql:before { content: 'SQL'; }
pre.src-sqlite:before { content: 'SQLite'; }
/* additional languages in org.el's org-babel-load-languages alist */
pre.src-forth:before { content: 'Forth'; }
pre.src-io:before { content: 'IO'; }
pre.src-J:before { content: 'J'; }
pre.src-makefile:before { content: 'Makefile'; }
pre.src-maxima:before { content: 'Maxima'; }
pre.src-perl:before { content: 'Perl'; }
pre.src-picolisp:before { content: 'Pico Lisp'; }
pre.src-scala:before { content: 'Scala'; }
pre.src-shell:before { content: 'Shell Script'; }
pre.src-ebnf2ps:before { content: 'ebfn2ps'; }
/* additional language identifiers per "defun org-babel-execute"
in ob-*.el */
pre.src-cpp:before { content: 'C++'; }
pre.src-abc:before { content: 'ABC'; }
pre.src-coq:before { content: 'Coq'; }
pre.src-groovy:before { content: 'Groovy'; }
/* additional language identifiers from org-babel-shell-names in
ob-shell.el: ob-shell is the only babel language using a lambda to put
the execution function name together. */
pre.src-bash:before { content: 'bash'; }
pre.src-csh:before { content: 'csh'; }
pre.src-ash:before { content: 'ash'; }
pre.src-dash:before { content: 'dash'; }
pre.src-ksh:before { content: 'ksh'; }
pre.src-mksh:before { content: 'mksh'; }
pre.src-posh:before { content: 'posh'; }
/* Additional Emacs modes also supported by the LaTeX listings package */
pre.src-ada:before { content: 'Ada'; }
pre.src-asm:before { content: 'Assembler'; }
pre.src-caml:before { content: 'Caml'; }
pre.src-delphi:before { content: 'Delphi'; }
pre.src-html:before { content: 'HTML'; }
pre.src-idl:before { content: 'IDL'; }
pre.src-mercury:before { content: 'Mercury'; }
pre.src-metapost:before { content: 'MetaPost'; }
pre.src-modula-2:before { content: 'Modula-2'; }
pre.src-pascal:before { content: 'Pascal'; }
pre.src-ps:before { content: 'PostScript'; }
pre.src-prolog:before { content: 'Prolog'; }
pre.src-simula:before { content: 'Simula'; }
pre.src-tcl:before { content: 'tcl'; }
pre.src-tex:before { content: 'TeX'; }
pre.src-plain-tex:before { content: 'Plain TeX'; }
pre.src-verilog:before { content: 'Verilog'; }
pre.src-vhdl:before { content: 'VHDL'; }
pre.src-xml:before { content: 'XML'; }
pre.src-nxml:before { content: 'XML'; }
/* add a generic configuration mode; LaTeX export needs an additional
(add-to-list 'org-latex-listings-langs '(conf " ")) in .emacs */
pre.src-conf:before { content: 'Configuration File'; }
table { border-collapse:collapse; }
caption.t-above { caption-side: top; }
caption.t-bottom { caption-side: bottom; }
td, th { vertical-align:top; }
th.org-right { text-align: center; }
th.org-left { text-align: center; }
th.org-center { text-align: center; }
td.org-right { text-align: right; }
td.org-left { text-align: left; }
td.org-center { text-align: center; }
dt { font-weight: bold; }
.footpara { display: inline; }
.footdef { margin-bottom: 1em; }
.figure { padding: 1em; }
.figure p { text-align: center; }
.inlinetask {
padding: 10px;
border: 2px solid gray;
margin: 10px;
background: #ffffcc;
}
#org-div-home-and-up
{ text-align: right; font-size: 70%; white-space: nowrap; }
textarea { overflow-x: auto; }
.linenr { font-size: smaller }
.code-highlighted { background-color: #ffff00; }
.org-info-js_info-navigation { border-style: none; }
#org-info-js_console-label
{ font-size: 10px; font-weight: bold; white-space: nowrap; }
.org-info-js_search-highlight
{ background-color: #ffff00; color: #000000; font-weight: bold; }
.org-svg { width: 90%; }
/*]]>*/-->
</style>
<link rel="stylesheet" title="Standard" href="qmckl.css" type="text/css" />
<script type="text/javascript" src="org-info.js">
/**
*
* @source: org-info.js
*
* @licstart The following is the entire license notice for the
* JavaScript code in org-info.js.
*
* Copyright (C) 2012-2019 Free Software Foundation, Inc.
*
*
* The JavaScript code in this tag is free software: you can
* redistribute it and/or modify it under the terms of the GNU
* General Public License (GNU GPL) as published by the Free Software
* Foundation, either version 3 of the License, or (at your option)
* any later version. The code is distributed WITHOUT ANY WARRANTY;
* without even the implied warranty of MERCHANTABILITY or FITNESS
* FOR A PARTICULAR PURPOSE. See the GNU GPL for more details.
*
* As additional permission under GNU GPL version 3 section 7, you
* may distribute non-source (e.g., minimized or compacted) forms of
* that code without the copy of the GNU GPL normally required by
* section 4, provided you include this license notice and a URL
* through which recipients can access the Corresponding Source.
*
* @licend The above is the entire license notice
* for the JavaScript code in org-info.js.
*
*/
</script>
<script type="text/javascript">
/*
@licstart The following is the entire license notice for the
JavaScript code in this tag.
Copyright (C) 2012-2019 Free Software Foundation, Inc.
The JavaScript code in this tag is free software: you can
redistribute it and/or modify it under the terms of the GNU
General Public License (GNU GPL) as published by the Free Software
Foundation, either version 3 of the License, or (at your option)
any later version. The code is distributed WITHOUT ANY WARRANTY;
without even the implied warranty of MERCHANTABILITY or FITNESS
FOR A PARTICULAR PURPOSE. See the GNU GPL for more details.
As additional permission under GNU GPL version 3 section 7, you
may distribute non-source (e.g., minimized or compacted) forms of
that code without the copy of the GNU GPL normally required by
section 4, provided you include this license notice and a URL
through which recipients can access the Corresponding Source.
@licend The above is the entire license notice
for the JavaScript code in this tag.
*/
<!--/*--><![CDATA[/*><!--*/
org_html_manager.set("TOC_DEPTH", "4");
org_html_manager.set("LINK_HOME", "index.html");
org_html_manager.set("LINK_UP", "");
org_html_manager.set("LOCAL_TOC", "1");
org_html_manager.set("VIEW_BUTTONS", "0");
org_html_manager.set("MOUSE_HINT", "underline");
org_html_manager.set("FIXED_TOC", "0");
org_html_manager.set("TOC", "1");
org_html_manager.set("VIEW", "info");
org_html_manager.setup(); // activate after the parameters are set
/*]]>*///-->
</script>
<script type="text/javascript">
/*
@licstart The following is the entire license notice for the
JavaScript code in this tag.
Copyright (C) 2012-2019 Free Software Foundation, Inc.
The JavaScript code in this tag is free software: you can
redistribute it and/or modify it under the terms of the GNU
General Public License (GNU GPL) as published by the Free Software
Foundation, either version 3 of the License, or (at your option)
any later version. The code is distributed WITHOUT ANY WARRANTY;
without even the implied warranty of MERCHANTABILITY or FITNESS
FOR A PARTICULAR PURPOSE. See the GNU GPL for more details.
As additional permission under GNU GPL version 3 section 7, you
may distribute non-source (e.g., minimized or compacted) forms of
that code without the copy of the GNU GPL normally required by
section 4, provided you include this license notice and a URL
through which recipients can access the Corresponding Source.
@licend The above is the entire license notice
for the JavaScript code in this tag.
*/
<!--/*--><![CDATA[/*><!--*/
function CodeHighlightOn(elem, id)
{
var target = document.getElementById(id);
if(null != target) {
elem.cacheClassElem = elem.className;
elem.cacheClassTarget = target.className;
target.className = "code-highlighted";
elem.className = "code-highlighted";
}
}
function CodeHighlightOff(elem, id)
{
var target = document.getElementById(id);
if(elem.cacheClassElem)
elem.className = elem.cacheClassElem;
if(elem.cacheClassTarget)
target.className = elem.cacheClassTarget;
}
/*]]>*///-->
</script>
<script type="text/x-mathjax-config">
MathJax.Hub.Config({
displayAlign: "center",
displayIndent: "0em",
"HTML-CSS": { scale: 100,
linebreaks: { automatic: "false" },
webFont: "TeX"
},
SVG: {scale: 100,
linebreaks: { automatic: "false" },
font: "TeX"},
NativeMML: {scale: 100},
TeX: { equationNumbers: {autoNumber: "AMS"},
MultLineWidth: "85%",
TagSide: "right",
TagIndent: ".8em"
}
});
</script>
<script type="text/javascript"
src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.0/MathJax.js?config=TeX-AMS_HTML"></script>
</head>
<body>
<div id="org-div-home-and-up">
<a accesskey="h" href=""> UP </a>
|
<a accesskey="H" href="index.html"> HOME </a>
</div><div id="content">
<h1 class="title">BLAS functions</h1>
<div id="table-of-contents">
<h2>Table of Contents</h2>
<div id="text-table-of-contents">
<ul>
<li><a href="#org2441856">1. Matrix operations</a>
<ul>
<li><a href="#orgb434b3c">1.1. <code>qmckl_dgemm</code></a></li>
<li><a href="#orgbd53465">1.2. <code>qmckl_adjugate</code></a></li>
</ul>
</li>
</ul>
</div>
</div>
<div id="outline-container-org2441856" class="outline-2">
<h2 id="org2441856"><span class="section-number-2">1</span> Matrix operations</h2>
<div class="outline-text-2" id="text-1">
</div>
<div id="outline-container-orgb434b3c" class="outline-3">
<h3 id="orgb434b3c"><span class="section-number-3">1.1</span> <code>qmckl_dgemm</code></h3>
<div class="outline-text-3" id="text-1-1">
<p>
Matrix multiplication:
</p>
<p>
\[
C_{ij} = \beta C_{ij} + \alpha \sum_{k} A_{ik} \cdot B_{kj}
\]
</p>
<table id="org798334b" border="2" cellspacing="0" cellpadding="6" rules="groups" frame="hsides">
<colgroup>
<col class="org-left" />
<col class="org-left" />
<col class="org-left" />
<col class="org-left" />
</colgroup>
<thead>
<tr>
<th scope="col" class="org-left">Variable</th>
<th scope="col" class="org-left">Type</th>
<th scope="col" class="org-left">In/Out</th>
<th scope="col" class="org-left">Description</th>
</tr>
</thead>
<tbody>
<tr>
<td class="org-left"><code>context</code></td>
<td class="org-left"><code>qmckl_context</code></td>
<td class="org-left">in</td>
<td class="org-left">Global state</td>
</tr>
<tr>
<td class="org-left"><code>TransA</code></td>
<td class="org-left"><code>char</code></td>
<td class="org-left">in</td>
<td class="org-left">'T' is transposed</td>
</tr>
<tr>
<td class="org-left"><code>TransB</code></td>
<td class="org-left"><code>char</code></td>
<td class="org-left">in</td>
<td class="org-left">'T' is transposed</td>
</tr>
<tr>
<td class="org-left"><code>m</code></td>
<td class="org-left"><code>int64_t</code></td>
<td class="org-left">in</td>
<td class="org-left">Number of rows of the input matrix</td>
</tr>
<tr>
<td class="org-left"><code>n</code></td>
<td class="org-left"><code>int64_t</code></td>
<td class="org-left">in</td>
<td class="org-left">Number of columns of the input matrix</td>
</tr>
<tr>
<td class="org-left"><code>k</code></td>
<td class="org-left"><code>int64_t</code></td>
<td class="org-left">in</td>
<td class="org-left">Number of columns of the input matrix</td>
</tr>
<tr>
<td class="org-left"><code>alpha</code></td>
<td class="org-left"><code>double</code></td>
<td class="org-left">in</td>
<td class="org-left">Number of columns of the input matrix</td>
</tr>
<tr>
<td class="org-left"><code>A</code></td>
<td class="org-left"><code>double[][lda]</code></td>
<td class="org-left">in</td>
<td class="org-left">Array containing the matrix \(A\)</td>
</tr>
<tr>
<td class="org-left"><code>lda</code></td>
<td class="org-left"><code>int64_t</code></td>
<td class="org-left">in</td>
<td class="org-left">Leading dimension of array <code>A</code></td>
</tr>
<tr>
<td class="org-left"><code>B</code></td>
<td class="org-left"><code>double[][ldb]</code></td>
<td class="org-left">in</td>
<td class="org-left">Array containing the matrix \(B\)</td>
</tr>
<tr>
<td class="org-left"><code>ldb</code></td>
<td class="org-left"><code>int64_t</code></td>
<td class="org-left">in</td>
<td class="org-left">Leading dimension of array <code>B</code></td>
</tr>
<tr>
<td class="org-left"><code>beta</code></td>
<td class="org-left"><code>double</code></td>
<td class="org-left">in</td>
<td class="org-left">Array containing the matrix \(B\)</td>
</tr>
<tr>
<td class="org-left"><code>C</code></td>
<td class="org-left"><code>double[][ldc]</code></td>
<td class="org-left">out</td>
<td class="org-left">Array containing the matrix \(B\)</td>
</tr>
<tr>
<td class="org-left"><code>ldc</code></td>
<td class="org-left"><code>int64_t</code></td>
<td class="org-left">in</td>
<td class="org-left">Leading dimension of array <code>B</code></td>
</tr>
</tbody>
</table>
<p>
Requirements:
</p>
<ul class="org-ul">
<li><code>context</code> is not <code>QMCKL_NULL_CONTEXT</code></li>
<li><code>m &gt; 0</code></li>
<li><code>n &gt; 0</code></li>
<li><code>k &gt; 0</code></li>
<li><code>lda &gt;= m</code></li>
<li><code>ldb &gt;= n</code></li>
<li><code>ldc &gt;= n</code></li>
<li><code>A</code> is allocated with at least \(m \times k \times 8\) bytes</li>
<li><code>B</code> is allocated with at least \(k \times n \times 8\) bytes</li>
<li><code>C</code> is allocated with at least \(m \times n \times 8\) bytes</li>
</ul>
<div class="org-src-container">
<pre class="src src-c"><span style="color: #228b22;">qmckl_exit_code</span> <span style="color: #0000ff;">qmckl_dgemm</span> (
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">qmckl_context</span> <span style="color: #a0522d;">context</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">char</span> <span style="color: #a0522d;">TransA</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">char</span> <span style="color: #a0522d;">TransB</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">int64_t</span> <span style="color: #a0522d;">m</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">int64_t</span> <span style="color: #a0522d;">n</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">int64_t</span> <span style="color: #a0522d;">k</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">double</span> <span style="color: #a0522d;">alpha</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">double</span>* <span style="color: #a0522d;">A</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">int64_t</span> <span style="color: #a0522d;">lda</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">double</span>* <span style="color: #a0522d;">B</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">int64_t</span> <span style="color: #a0522d;">ldb</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">double</span> <span style="color: #a0522d;">beta</span>,
<span style="color: #228b22;">double</span>* <span style="color: #a020f0;">const</span> <span style="color: #a0522d;">C</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">int64_t</span> <span style="color: #a0522d;">ldc</span> );
</pre>
</div>
<div class="org-src-container">
<pre class="src src-f90"><span style="color: #228b22;">integer</span><span style="color: #a0522d;"> function qmckl_dgemm_f(context, TransA, TransB, </span><span style="color: #a020f0;">&amp;</span>
m, n, k, alpha, A, LDA, B, LDB, beta, C, LDC) <span style="color: #a020f0;">&amp;</span>
<span style="color: #a020f0;">result</span>(info)
<span style="color: #a020f0;">use</span> <span style="color: #0000ff;">qmckl</span>
<span style="color: #a020f0;">implicit</span> <span style="color: #228b22;">none</span>
<span style="color: #228b22;">integer</span>(qmckl_context), <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> context</span>
<span style="color: #228b22;">character</span> , <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> TransA, TransB</span>
<span style="color: #228b22;">integer</span>*8 , <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> m, n, k</span>
<span style="color: #228b22;">double precision</span> , <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> alpha, beta</span>
<span style="color: #228b22;">integer</span>*8 , <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> lda</span>
<span style="color: #228b22;">double precision</span> , <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> A(lda,*)</span>
<span style="color: #228b22;">integer</span>*8 , <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> ldb</span>
<span style="color: #228b22;">double precision</span> , <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> B(ldb,*)</span>
<span style="color: #228b22;">integer</span>*8 , <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> ldc</span>
<span style="color: #228b22;">double precision</span> , <span style="color: #a020f0;">intent</span>(out) ::<span style="color: #a0522d;"> C(ldc,*)</span>
info = QMCKL_SUCCESS
<span style="color: #a020f0;">if</span> (context == QMCKL_NULL_CONTEXT) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_CONTEXT
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">if</span> (m &lt;= 0_8) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_ARG_4
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">if</span> (n &lt;= 0_8) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_ARG_5
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">if</span> (k &lt;= 0_8) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_ARG_6
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">if</span> (LDA &lt;= 0) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_ARG_9
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">if</span> (LDB &lt;= 0) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_ARG_11
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">if</span> (LDC &lt;= 0) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_ARG_14
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">call</span> <span style="color: #0000ff;">dgemm</span>(transA, transB, <span style="color: #a020f0;">int</span>(m,4), <span style="color: #a020f0;">int</span>(n,4), <span style="color: #a020f0;">int</span>(k,4), <span style="color: #a020f0;">&amp;</span>
alpha, A, <span style="color: #a020f0;">int</span>(LDA,4), B, <span style="color: #a020f0;">int</span>(LDB,4), beta, C, <span style="color: #a020f0;">int</span>(LDC,4))
<span style="color: #a020f0;">end function</span> <span style="color: #0000ff;">qmckl_dgemm_f</span>
</pre>
</div>
</div>
</div>
<div id="outline-container-orgbd53465" class="outline-3">
<h3 id="orgbd53465"><span class="section-number-3">1.2</span> <code>qmckl_adjugate</code></h3>
<div class="outline-text-3" id="text-1-2">
<p>
Given a matrix \(\mathbf{A}\), the adjugate matrix
\(\text{adj}(\mathbf{A})\) is the transpose of the cofactors matrix
of \(\mathbf{A}\).
</p>
<p>
\[
\mathbf{B} = \text{adj}(\mathbf{A}) = \text{det}(\mathbf{A}) \, \mathbf{A}^{-1}
\]
</p>
<p>
See also: <a href="https://en.wikipedia.org/wiki/Adjugate_matrix">https://en.wikipedia.org/wiki/Adjugate_matrix</a>
</p>
<table id="org82ea90c" border="2" cellspacing="0" cellpadding="6" rules="groups" frame="hsides">
<colgroup>
<col class="org-left" />
<col class="org-left" />
<col class="org-left" />
<col class="org-left" />
</colgroup>
<thead>
<tr>
<th scope="col" class="org-left">Variable</th>
<th scope="col" class="org-left">Type</th>
<th scope="col" class="org-left">In/Out</th>
<th scope="col" class="org-left">Description</th>
</tr>
</thead>
<tbody>
<tr>
<td class="org-left"><code>context</code></td>
<td class="org-left"><code>qmckl_context</code></td>
<td class="org-left">in</td>
<td class="org-left">Global state</td>
</tr>
<tr>
<td class="org-left"><code>n</code></td>
<td class="org-left"><code>int64_t</code></td>
<td class="org-left">in</td>
<td class="org-left">Number of rows and columns of the input matrix</td>
</tr>
<tr>
<td class="org-left"><code>A</code></td>
<td class="org-left"><code>double[][lda]</code></td>
<td class="org-left">in</td>
<td class="org-left">Array containing the \(n \times n\) matrix \(A\)</td>
</tr>
<tr>
<td class="org-left"><code>lda</code></td>
<td class="org-left"><code>int64_t</code></td>
<td class="org-left">in</td>
<td class="org-left">Leading dimension of array <code>A</code></td>
</tr>
<tr>
<td class="org-left"><code>B</code></td>
<td class="org-left"><code>double[][ldb]</code></td>
<td class="org-left">out</td>
<td class="org-left">Adjugate of \(A\)</td>
</tr>
<tr>
<td class="org-left"><code>ldb</code></td>
<td class="org-left"><code>int64_t</code></td>
<td class="org-left">in</td>
<td class="org-left">Leading dimension of array <code>B</code></td>
</tr>
<tr>
<td class="org-left"><code>det_l</code></td>
<td class="org-left"><code>double</code></td>
<td class="org-left">inout</td>
<td class="org-left">determinant of \(A\)</td>
</tr>
</tbody>
</table>
<p>
Requirements:
</p>
<ul class="org-ul">
<li><code>context</code> is not <code>QMCKL_NULL_CONTEXT</code></li>
<li><code>n &gt; 0</code></li>
<li><code>lda &gt;= m</code></li>
<li><code>A</code> is allocated with at least \(m \times m \times 8\) bytes</li>
<li><code>ldb &gt;= m</code></li>
<li><code>B</code> is allocated with at least \(m \times m \times 8\) bytes</li>
</ul>
<div class="org-src-container">
<pre class="src src-c"><span style="color: #228b22;">qmckl_exit_code</span> <span style="color: #0000ff;">qmckl_adjugate</span> (
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">qmckl_context</span> <span style="color: #a0522d;">context</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">int64_t</span> <span style="color: #a0522d;">n</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">double</span>* <span style="color: #a0522d;">A</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">int64_t</span> <span style="color: #a0522d;">lda</span>,
<span style="color: #228b22;">double</span>* <span style="color: #a020f0;">const</span> <span style="color: #a0522d;">B</span>,
<span style="color: #a020f0;">const</span> <span style="color: #228b22;">int64_t</span> <span style="color: #a0522d;">ldb</span>,
<span style="color: #228b22;">double</span>* <span style="color: #a0522d;">det_l</span> );
</pre>
</div>
<p>
For small matrices (&le; 5&times; 5), we use specialized routines
for performance motivations. For larger sizes, we rely on the
LAPACK library.
</p>
<div class="org-src-container">
<pre class="src src-f90"><span style="color: #228b22;">integer </span><span style="color: #a020f0;">function</span><span style="color: #a0522d;"> </span><span style="color: #0000ff;">qmckl_adjugate_f</span><span style="color: #000000; background-color: #ffffff;">(context, na, A, LDA, B, ldb, det_l)</span><span style="color: #a0522d;"> </span><span style="color: #a020f0;">&amp;</span>
<span style="color: #a020f0;">result</span>(info)
<span style="color: #a020f0;">use</span> <span style="color: #0000ff;">qmckl</span>
<span style="color: #a020f0;">implicit</span> <span style="color: #228b22;">none</span>
<span style="color: #228b22;">integer</span>(qmckl_context) , <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> context</span>
<span style="color: #228b22;">double precision</span>, <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> A (LDA,*)</span>
<span style="color: #228b22;">integer</span>*8, <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> LDA</span>
<span style="color: #228b22;">double precision</span>, <span style="color: #a020f0;">intent</span>(out) ::<span style="color: #a0522d;"> B (LDB,*)</span>
<span style="color: #228b22;">integer</span>*8, <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> LDB</span>
<span style="color: #228b22;">integer</span>*8, <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> na</span>
<span style="color: #228b22;">double precision</span>, <span style="color: #a020f0;">intent</span>(inout) ::<span style="color: #a0522d;"> det_l</span>
<span style="color: #228b22;">integer</span> ::<span style="color: #a0522d;"> i,j</span>
info = QMCKL_SUCCESS
<span style="color: #a020f0;">if</span> (context == QMCKL_NULL_CONTEXT) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_CONTEXT
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">if</span> (na &lt;= 0_8) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_ARG_2
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">if</span> (LDA &lt;= 0_8) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_ARG_4
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">if</span> (LDA &lt; na) <span style="color: #a020f0;">then</span>
info = QMCKL_INVALID_ARG_4
<span style="color: #a020f0;">return</span>
<span style="color: #a020f0;">endif</span>
<span style="color: #a020f0;">select case</span> (na)
<span style="color: #a020f0;">case</span> (5)
<span style="color: #a020f0;">call</span> <span style="color: #0000ff;">adjugate5</span>(A,LDA,B,LDB,na,det_l)
<span style="color: #a020f0;">case</span> (4)
<span style="color: #a020f0;">call</span> <span style="color: #0000ff;">adjugate4</span>(A,LDA,B,LDB,na,det_l)
<span style="color: #a020f0;">case</span> (3)
<span style="color: #a020f0;">call</span> <span style="color: #0000ff;">adjugate3</span>(A,LDA,B,LDB,na,det_l)
<span style="color: #a020f0;">case</span> (2)
<span style="color: #a020f0;">call</span> <span style="color: #0000ff;">adjugate2</span>(A,LDA,B,LDB,na,det_l)
<span style="color: #a020f0;">case</span> (1)
det_l = a(1,1)
b(1,1) = 1.d0
<span style="color: #a020f0;">case</span> default
<span style="color: #a020f0;">call</span> <span style="color: #0000ff;">adjugate_general</span>(context, na, A, LDA, B, LDB, det_l)
<span style="color: #a020f0;">end select</span>
<span style="color: #a020f0;">end function</span> <span style="color: #0000ff;">qmckl_adjugate_f</span>
</pre>
</div>
<div class="org-src-container">
<pre class="src src-f90"><span style="color: #a020f0;">subroutine</span> <span style="color: #0000ff;">adjugate_general</span>(context, na, A, LDA, B, LDB, det_l)
<span style="color: #a020f0;">use</span> <span style="color: #0000ff;">qmckl</span>
<span style="color: #a020f0;">implicit</span> <span style="color: #228b22;">none</span>
<span style="color: #228b22;">integer</span>(qmckl_context), <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> context</span>
<span style="color: #228b22;">double precision</span>, <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> A (LDA,na)</span>
<span style="color: #228b22;">integer</span>*8, <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> LDA</span>
<span style="color: #228b22;">double precision</span>, <span style="color: #a020f0;">intent</span>(out) ::<span style="color: #a0522d;"> B (LDB,na)</span>
<span style="color: #228b22;">integer</span>*8, <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> LDB</span>
<span style="color: #228b22;">integer</span>*8, <span style="color: #a020f0;">intent</span>(in) ::<span style="color: #a0522d;"> na</span>
<span style="color: #228b22;">double precision</span>, <span style="color: #a020f0;">intent</span>(inout) ::<span style="color: #a0522d;"> det_l</span>
<span style="color: #228b22;">double precision</span> ::<span style="color: #a0522d;"> work(LDA*max(na,64))</span>
<span style="color: #228b22;">integer</span> ::<span style="color: #a0522d;"> inf</span>
<span style="color: #228b22;">integer</span> ::<span style="color: #a0522d;"> ipiv(LDA)</span>
<span style="color: #228b22;">integer</span> ::<span style="color: #a0522d;"> lwork</span>
<span style="color: #228b22;">integer</span> ::<span style="color: #a0522d;"> i, j</span>
</pre>
</div>
<p>
We first copy the array <code>A</code> into array <code>B</code>.
</p>
<div class="org-src-container">
<pre class="src src-f90">B(1:na,1:na) = A(1:na,1:na)
</pre>
</div>
<p>
Then, we compute the LU factorization \(LU=A\)
using the <code>dgetrf</code> routine.
</p>
<div class="org-src-container">
<pre class="src src-f90"><span style="color: #a020f0;">call</span> <span style="color: #0000ff;">dgetrf</span>(na, na, B, LDB, ipiv, inf )
</pre>
</div>
<p>
By convention, the determinant of \(L\) is equal to one, so the
determinant of \(A\) is equal to the determinant of \(U\), which is
simply computed as the product of its diagonal elements.
</p>
<div class="org-src-container">
<pre class="src src-f90">det_l = 1.d0
j=0
<span style="color: #a020f0;">do</span> i=1,na
j = j+<span style="color: #a020f0;">min</span>(<span style="color: #a020f0;">abs</span>(ipiv(i)-i),1)
det_l = det_l*B(i,i)
<span style="color: #a020f0;">enddo</span>
</pre>
</div>
<p>
As <code>dgetrf</code> returns \(PLU=A\) where \(P\) is a permutation matrix, the
sign of the determinant is computed as \(-1^m\) where \(m\) is the
number of permutations.
</p>
<div class="org-src-container">
<pre class="src src-f90"><span style="color: #a020f0;">if</span> (<span style="color: #a020f0;">iand</span>(j,1) /= 0) <span style="color: #a020f0;">then</span>
det_l = -det_l
<span style="color: #a020f0;">endif</span>
</pre>
</div>
<p>
Then, the inverse of \(A\) is computed using <code>dgetri</code>:
</p>
<div class="org-src-container">
<pre class="src src-f90">lwork = <span style="color: #a020f0;">SIZE</span>(work)
<span style="color: #a020f0;">call</span> <span style="color: #0000ff;">dgetri</span>(na, B, LDB, ipiv, work, lwork, inf )
</pre>
</div>
<p>
and the adjugate matrix is computed as the product of the
determinant with the inverse:
</p>
<div class="org-src-container">
<pre class="src src-f90"> B(:,:) = B(:,:)*det_l
<span style="color: #a020f0;">end subroutine</span> <span style="color: #0000ff;">adjugate_general</span>
</pre>
</div>
</div>
</div>
</div>
</div>
<div id="postamble" class="status">
<p class="author">Author: TREX CoE</p>
<p class="date">Created: 2022-01-08 Sat 14:20</p>
<p class="validation"><a href="http://validator.w3.org/check?uri=referer">Validate</a></p>
</div>
</body>
</html>