From 511f1e51c009a51be8a2839a8b0bb0e81106bdf3 Mon Sep 17 00:00:00 2001 From: caremoli Date: Wed, 26 Sep 2007 14:45:33 +0000 Subject: [PATCH 1/1] initial import into CVS --- AUTHORS | 0 COPYING | 340 +++ ChangeLog | 0 INSTALL | 236 ++ Makefile.am | 10 + NEWS | 0 README | 19 + adm_local/check_Kernel.m4 | 59 + adm_local/check_omniorb.m4 | 251 ++ adm_local/make_common_starter.am | 11 + adm_local/python.m4 | 163 ++ autogen.sh | 12 + configure.ac | 37 + idl/DSCCODE.idl | 55 + idl/Makefile.am | 69 + resources/DSCCODECatalog.xml | 330 +++ resources/Makefile.am | 9 + resources/calcium4.xml | 209 ++ resources/stream2.xml | 51 + src/DSCCODAENG/DSCCODAENG.cxx | 141 ++ src/DSCCODAENG/DSCCODAENG.hxx | 34 + src/DSCCODAENG/Makefile.am | 9 + src/DSCCODBENG/DSCCODBENG.cxx | 144 ++ src/DSCCODBENG/DSCCODBENG.hxx | 34 + src/DSCCODBENG/Makefile.am | 9 + src/DSCCODCENG/DSCCODCENG.cxx | 162 ++ src/DSCCODCENG/DSCCODCENG.hxx | 29 + src/DSCCODCENG/Makefile.am | 9 + src/DSCCODDENG/DSCCODDENG.cxx | 160 ++ src/DSCCODDENG/DSCCODDENG.hxx | 31 + src/DSCCODDENG/Makefile.am | 9 + src/FLUIDE/FLUIDE.cxx | 154 ++ src/FLUIDE/FLUIDE.hxx | 29 + src/FLUIDE/Makefile.am | 12 + src/FLUIDE/fluide.f | 161 ++ src/INTERPI/INTERPI.cxx | 175 ++ src/INTERPI/INTERPI.hxx | 27 + src/INTERPI/Makefile.am | 15 + src/INTERPI/inter.f | 45 + src/INTERPI/interi.f | 39 + src/Makefile.am | 1 + src/NEUTRO/Makefile.am | 12 + src/NEUTRO/NEUTRO.cxx | 157 ++ src/NEUTRO/NEUTRO.hxx | 29 + src/NEUTRO/neutro.f | 162 ++ src/SOLIDE/Makefile.am | 12 + src/SOLIDE/SOLIDE.cxx | 164 ++ src/SOLIDE/SOLIDE.hxx | 29 + src/SOLIDE/solid.f | 3801 ++++++++++++++++++++++++++++++ 49 files changed, 7656 insertions(+) create mode 100644 AUTHORS create mode 100644 COPYING create mode 100644 ChangeLog create mode 100644 INSTALL create mode 100644 Makefile.am create mode 100644 NEWS create mode 100644 README create mode 100644 adm_local/check_Kernel.m4 create mode 100644 adm_local/check_omniorb.m4 create mode 100644 adm_local/make_common_starter.am create mode 100644 adm_local/python.m4 create mode 100755 autogen.sh create mode 100644 configure.ac create mode 100644 idl/DSCCODE.idl create mode 100644 idl/Makefile.am create mode 100644 resources/DSCCODECatalog.xml create mode 100644 resources/Makefile.am create mode 100644 resources/calcium4.xml create mode 100644 resources/stream2.xml create mode 100755 src/DSCCODAENG/DSCCODAENG.cxx create mode 100644 src/DSCCODAENG/DSCCODAENG.hxx create mode 100644 src/DSCCODAENG/Makefile.am create mode 100755 src/DSCCODBENG/DSCCODBENG.cxx create mode 100644 src/DSCCODBENG/DSCCODBENG.hxx create mode 100644 src/DSCCODBENG/Makefile.am create mode 100755 src/DSCCODCENG/DSCCODCENG.cxx create mode 100644 src/DSCCODCENG/DSCCODCENG.hxx create mode 100644 src/DSCCODCENG/Makefile.am create mode 100755 src/DSCCODDENG/DSCCODDENG.cxx create mode 100644 src/DSCCODDENG/DSCCODDENG.hxx create mode 100644 src/DSCCODDENG/Makefile.am create mode 100755 src/FLUIDE/FLUIDE.cxx create mode 100644 src/FLUIDE/FLUIDE.hxx create mode 100644 src/FLUIDE/Makefile.am create mode 100644 src/FLUIDE/fluide.f create mode 100755 src/INTERPI/INTERPI.cxx create mode 100644 src/INTERPI/INTERPI.hxx create mode 100644 src/INTERPI/Makefile.am create mode 100644 src/INTERPI/inter.f create mode 100644 src/INTERPI/interi.f create mode 100644 src/Makefile.am create mode 100644 src/NEUTRO/Makefile.am create mode 100755 src/NEUTRO/NEUTRO.cxx create mode 100644 src/NEUTRO/NEUTRO.hxx create mode 100644 src/NEUTRO/neutro.f create mode 100644 src/SOLIDE/Makefile.am create mode 100755 src/SOLIDE/SOLIDE.cxx create mode 100644 src/SOLIDE/SOLIDE.hxx create mode 100644 src/SOLIDE/solid.f diff --git a/AUTHORS b/AUTHORS new file mode 100644 index 0000000..e69de29 diff --git a/COPYING b/COPYING new file mode 100644 index 0000000..d60c31a --- /dev/null +++ b/COPYING @@ -0,0 +1,340 @@ + GNU GENERAL PUBLIC LICENSE + Version 2, June 1991 + + Copyright (C) 1989, 1991 Free Software Foundation, Inc. + 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + Everyone is permitted to copy and distribute verbatim copies + of this license document, but changing it is not allowed. + + Preamble + + The licenses for most software are designed to take away your +freedom to share and change it. By contrast, the GNU General Public +License is intended to guarantee your freedom to share and change free +software--to make sure the software is free for all its users. This +General Public License applies to most of the Free Software +Foundation's software and to any other program whose authors commit to +using it. (Some other Free Software Foundation software is covered by +the GNU Library General Public License instead.) You can apply it to +your programs, too. + + When we speak of free software, we are referring to freedom, not +price. Our General Public Licenses are designed to make sure that you +have the freedom to distribute copies of free software (and charge for +this service if you wish), that you receive source code or can get it +if you want it, that you can change the software or use pieces of it +in new free programs; and that you know you can do these things. + + To protect your rights, we need to make restrictions that forbid +anyone to deny you these rights or to ask you to surrender the rights. +These restrictions translate to certain responsibilities for you if you +distribute copies of the software, or if you modify it. + + For example, if you distribute copies of such a program, whether +gratis or for a fee, you must give the recipients all the rights that +you have. You must make sure that they, too, receive or can get the +source code. And you must show them these terms so they know their +rights. + + We protect your rights with two steps: (1) copyright the software, and +(2) offer you this license which gives you legal permission to copy, +distribute and/or modify the software. + + Also, for each author's protection and ours, we want to make certain +that everyone understands that there is no warranty for this free +software. If the software is modified by someone else and passed on, we +want its recipients to know that what they have is not the original, so +that any problems introduced by others will not reflect on the original +authors' reputations. + + Finally, any free program is threatened constantly by software +patents. We wish to avoid the danger that redistributors of a free +program will individually obtain patent licenses, in effect making the +program proprietary. To prevent this, we have made it clear that any +patent must be licensed for everyone's free use or not licensed at all. + + The precise terms and conditions for copying, distribution and +modification follow. + + GNU GENERAL PUBLIC LICENSE + TERMS AND CONDITIONS FOR COPYING, DISTRIBUTION AND MODIFICATION + + 0. This License applies to any program or other work which contains +a notice placed by the copyright holder saying it may be distributed +under the terms of this General Public License. The "Program", below, +refers to any such program or work, and a "work based on the Program" +means either the Program or any derivative work under copyright law: +that is to say, a work containing the Program or a portion of it, +either verbatim or with modifications and/or translated into another +language. (Hereinafter, translation is included without limitation in +the term "modification".) Each licensee is addressed as "you". + +Activities other than copying, distribution and modification are not +covered by this License; they are outside its scope. The act of +running the Program is not restricted, and the output from the Program +is covered only if its contents constitute a work based on the +Program (independent of having been made by running the Program). +Whether that is true depends on what the Program does. + + 1. You may copy and distribute verbatim copies of the Program's +source code as you receive it, in any medium, provided that you +conspicuously and appropriately publish on each copy an appropriate +copyright notice and disclaimer of warranty; keep intact all the +notices that refer to this License and to the absence of any warranty; +and give any other recipients of the Program a copy of this License +along with the Program. + +You may charge a fee for the physical act of transferring a copy, and +you may at your option offer warranty protection in exchange for a fee. + + 2. You may modify your copy or copies of the Program or any portion +of it, thus forming a work based on the Program, and copy and +distribute such modifications or work under the terms of Section 1 +above, provided that you also meet all of these conditions: + + a) You must cause the modified files to carry prominent notices + stating that you changed the files and the date of any change. + + b) You must cause any work that you distribute or publish, that in + whole or in part contains or is derived from the Program or any + part thereof, to be licensed as a whole at no charge to all third + parties under the terms of this License. + + c) If the modified program normally reads commands interactively + when run, you must cause it, when started running for such + interactive use in the most ordinary way, to print or display an + announcement including an appropriate copyright notice and a + notice that there is no warranty (or else, saying that you provide + a warranty) and that users may redistribute the program under + these conditions, and telling the user how to view a copy of this + License. (Exception: if the Program itself is interactive but + does not normally print such an announcement, your work based on + the Program is not required to print an announcement.) + +These requirements apply to the modified work as a whole. If +identifiable sections of that work are not derived from the Program, +and can be reasonably considered independent and separate works in +themselves, then this License, and its terms, do not apply to those +sections when you distribute them as separate works. But when you +distribute the same sections as part of a whole which is a work based +on the Program, the distribution of the whole must be on the terms of +this License, whose permissions for other licensees extend to the +entire whole, and thus to each and every part regardless of who wrote it. + +Thus, it is not the intent of this section to claim rights or contest +your rights to work written entirely by you; rather, the intent is to +exercise the right to control the distribution of derivative or +collective works based on the Program. + +In addition, mere aggregation of another work not based on the Program +with the Program (or with a work based on the Program) on a volume of +a storage or distribution medium does not bring the other work under +the scope of this License. + + 3. You may copy and distribute the Program (or a work based on it, +under Section 2) in object code or executable form under the terms of +Sections 1 and 2 above provided that you also do one of the following: + + a) Accompany it with the complete corresponding machine-readable + source code, which must be distributed under the terms of Sections + 1 and 2 above on a medium customarily used for software interchange; or, + + b) Accompany it with a written offer, valid for at least three + years, to give any third party, for a charge no more than your + cost of physically performing source distribution, a complete + machine-readable copy of the corresponding source code, to be + distributed under the terms of Sections 1 and 2 above on a medium + customarily used for software interchange; or, + + c) Accompany it with the information you received as to the offer + to distribute corresponding source code. (This alternative is + allowed only for noncommercial distribution and only if you + received the program in object code or executable form with such + an offer, in accord with Subsection b above.) + +The source code for a work means the preferred form of the work for +making modifications to it. For an executable work, complete source +code means all the source code for all modules it contains, plus any +associated interface definition files, plus the scripts used to +control compilation and installation of the executable. However, as a +special exception, the source code distributed need not include +anything that is normally distributed (in either source or binary +form) with the major components (compiler, kernel, and so on) of the +operating system on which the executable runs, unless that component +itself accompanies the executable. + +If distribution of executable or object code is made by offering +access to copy from a designated place, then offering equivalent +access to copy the source code from the same place counts as +distribution of the source code, even though third parties are not +compelled to copy the source along with the object code. + + 4. You may not copy, modify, sublicense, or distribute the Program +except as expressly provided under this License. Any attempt +otherwise to copy, modify, sublicense or distribute the Program is +void, and will automatically terminate your rights under this License. +However, parties who have received copies, or rights, from you under +this License will not have their licenses terminated so long as such +parties remain in full compliance. + + 5. You are not required to accept this License, since you have not +signed it. However, nothing else grants you permission to modify or +distribute the Program or its derivative works. These actions are +prohibited by law if you do not accept this License. Therefore, by +modifying or distributing the Program (or any work based on the +Program), you indicate your acceptance of this License to do so, and +all its terms and conditions for copying, distributing or modifying +the Program or works based on it. + + 6. Each time you redistribute the Program (or any work based on the +Program), the recipient automatically receives a license from the +original licensor to copy, distribute or modify the Program subject to +these terms and conditions. You may not impose any further +restrictions on the recipients' exercise of the rights granted herein. +You are not responsible for enforcing compliance by third parties to +this License. + + 7. If, as a consequence of a court judgment or allegation of patent +infringement or for any other reason (not limited to patent issues), +conditions are imposed on you (whether by court order, agreement or +otherwise) that contradict the conditions of this License, they do not +excuse you from the conditions of this License. If you cannot +distribute so as to satisfy simultaneously your obligations under this +License and any other pertinent obligations, then as a consequence you +may not distribute the Program at all. For example, if a patent +license would not permit royalty-free redistribution of the Program by +all those who receive copies directly or indirectly through you, then +the only way you could satisfy both it and this License would be to +refrain entirely from distribution of the Program. + +If any portion of this section is held invalid or unenforceable under +any particular circumstance, the balance of the section is intended to +apply and the section as a whole is intended to apply in other +circumstances. + +It is not the purpose of this section to induce you to infringe any +patents or other property right claims or to contest validity of any +such claims; this section has the sole purpose of protecting the +integrity of the free software distribution system, which is +implemented by public license practices. Many people have made +generous contributions to the wide range of software distributed +through that system in reliance on consistent application of that +system; it is up to the author/donor to decide if he or she is willing +to distribute software through any other system and a licensee cannot +impose that choice. + +This section is intended to make thoroughly clear what is believed to +be a consequence of the rest of this License. + + 8. If the distribution and/or use of the Program is restricted in +certain countries either by patents or by copyrighted interfaces, the +original copyright holder who places the Program under this License +may add an explicit geographical distribution limitation excluding +those countries, so that distribution is permitted only in or among +countries not thus excluded. In such case, this License incorporates +the limitation as if written in the body of this License. + + 9. The Free Software Foundation may publish revised and/or new versions +of the General Public License from time to time. Such new versions will +be similar in spirit to the present version, but may differ in detail to +address new problems or concerns. + +Each version is given a distinguishing version number. If the Program +specifies a version number of this License which applies to it and "any +later version", you have the option of following the terms and conditions +either of that version or of any later version published by the Free +Software Foundation. If the Program does not specify a version number of +this License, you may choose any version ever published by the Free Software +Foundation. + + 10. If you wish to incorporate parts of the Program into other free +programs whose distribution conditions are different, write to the author +to ask for permission. For software which is copyrighted by the Free +Software Foundation, write to the Free Software Foundation; we sometimes +make exceptions for this. Our decision will be guided by the two goals +of preserving the free status of all derivatives of our free software and +of promoting the sharing and reuse of software generally. + + NO WARRANTY + + 11. BECAUSE THE PROGRAM IS LICENSED FREE OF CHARGE, THERE IS NO WARRANTY +FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE LAW. EXCEPT WHEN +OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR OTHER PARTIES +PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESSED +OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF +MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS +TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE +PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY SERVICING, +REPAIR OR CORRECTION. + + 12. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING +WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY MODIFY AND/OR +REDISTRIBUTE THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, +INCLUDING ANY GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING +OUT OF THE USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED +TO LOSS OF DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY +YOU OR THIRD PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER +PROGRAMS), EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE +POSSIBILITY OF SUCH DAMAGES. + + END OF TERMS AND CONDITIONS + + How to Apply These Terms to Your New Programs + + If you develop a new program, and you want it to be of the greatest +possible use to the public, the best way to achieve this is to make it +free software which everyone can redistribute and change under these terms. + + To do so, attach the following notices to the program. It is safest +to attach them to the start of each source file to most effectively +convey the exclusion of warranty; and each file should have at least +the "copyright" line and a pointer to where the full notice is found. + + + Copyright (C) + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program; if not, write to the Free Software + Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + + +Also add information on how to contact you by electronic and paper mail. + +If the program is interactive, make it output a short notice like this +when it starts in an interactive mode: + + Gnomovision version 69, Copyright (C) year name of author + Gnomovision comes with ABSOLUTELY NO WARRANTY; for details type `show w'. + This is free software, and you are welcome to redistribute it + under certain conditions; type `show c' for details. + +The hypothetical commands `show w' and `show c' should show the appropriate +parts of the General Public License. Of course, the commands you use may +be called something other than `show w' and `show c'; they could even be +mouse-clicks or menu items--whatever suits your program. + +You should also get your employer (if you work as a programmer) or your +school, if any, to sign a "copyright disclaimer" for the program, if +necessary. Here is a sample; alter the names: + + Yoyodyne, Inc., hereby disclaims all copyright interest in the program + `Gnomovision' (which makes passes at compilers) written by James Hacker. + + , 1 April 1989 + Ty Coon, President of Vice + +This General Public License does not permit incorporating your program into +proprietary programs. If your program is a subroutine library, you may +consider it more useful to permit linking proprietary applications with the +library. If this is what you want to do, use the GNU Library General +Public License instead of this License. diff --git a/ChangeLog b/ChangeLog new file mode 100644 index 0000000..e69de29 diff --git a/INSTALL b/INSTALL new file mode 100644 index 0000000..56b077d --- /dev/null +++ b/INSTALL @@ -0,0 +1,236 @@ +Installation Instructions +************************* + +Copyright (C) 1994, 1995, 1996, 1999, 2000, 2001, 2002, 2004, 2005 Free +Software Foundation, Inc. + +This file is free documentation; the Free Software Foundation gives +unlimited permission to copy, distribute and modify it. + +Basic Installation +================== + +These are generic installation instructions. + + The `configure' shell script attempts to guess correct values for +various system-dependent variables used during compilation. It uses +those values to create a `Makefile' in each directory of the package. +It may also create one or more `.h' files containing system-dependent +definitions. Finally, it creates a shell script `config.status' that +you can run in the future to recreate the current configuration, and a +file `config.log' containing compiler output (useful mainly for +debugging `configure'). + + It can also use an optional file (typically called `config.cache' +and enabled with `--cache-file=config.cache' or simply `-C') that saves +the results of its tests to speed up reconfiguring. (Caching is +disabled by default to prevent problems with accidental use of stale +cache files.) + + If you need to do unusual things to compile the package, please try +to figure out how `configure' could check whether to do them, and mail +diffs or instructions to the address given in the `README' so they can +be considered for the next release. If you are using the cache, and at +some point `config.cache' contains results you don't want to keep, you +may remove or edit it. + + The file `configure.ac' (or `configure.in') is used to create +`configure' by a program called `autoconf'. You only need +`configure.ac' if you want to change it or regenerate `configure' using +a newer version of `autoconf'. + +The simplest way to compile this package is: + + 1. `cd' to the directory containing the package's source code and type + `./configure' to configure the package for your system. If you're + using `csh' on an old version of System V, you might need to type + `sh ./configure' instead to prevent `csh' from trying to execute + `configure' itself. + + Running `configure' takes awhile. While running, it prints some + messages telling which features it is checking for. + + 2. Type `make' to compile the package. + + 3. Optionally, type `make check' to run any self-tests that come with + the package. + + 4. Type `make install' to install the programs and any data files and + documentation. + + 5. You can remove the program binaries and object files from the + source code directory by typing `make clean'. To also remove the + files that `configure' created (so you can compile the package for + a different kind of computer), type `make distclean'. There is + also a `make maintainer-clean' target, but that is intended mainly + for the package's developers. If you use it, you may have to get + all sorts of other programs in order to regenerate files that came + with the distribution. + +Compilers and Options +===================== + +Some systems require unusual options for compilation or linking that the +`configure' script does not know about. Run `./configure --help' for +details on some of the pertinent environment variables. + + You can give `configure' initial values for configuration parameters +by setting variables in the command line or in the environment. Here +is an example: + + ./configure CC=c89 CFLAGS=-O2 LIBS=-lposix + + *Note Defining Variables::, for more details. + +Compiling For Multiple Architectures +==================================== + +You can compile the package for more than one kind of computer at the +same time, by placing the object files for each architecture in their +own directory. To do this, you must use a version of `make' that +supports the `VPATH' variable, such as GNU `make'. `cd' to the +directory where you want the object files and executables to go and run +the `configure' script. `configure' automatically checks for the +source code in the directory that `configure' is in and in `..'. + + If you have to use a `make' that does not support the `VPATH' +variable, you have to compile the package for one architecture at a +time in the source code directory. After you have installed the +package for one architecture, use `make distclean' before reconfiguring +for another architecture. + +Installation Names +================== + +By default, `make install' will install the package's files in +`/usr/local/bin', `/usr/local/man', etc. You can specify an +installation prefix other than `/usr/local' by giving `configure' the +option `--prefix=PREFIX'. + + You can specify separate installation prefixes for +architecture-specific files and architecture-independent files. If you +give `configure' the option `--exec-prefix=PREFIX', the package will +use PREFIX as the prefix for installing programs and libraries. +Documentation and other data files will still use the regular prefix. + + In addition, if you use an unusual directory layout you can give +options like `--bindir=DIR' to specify different values for particular +kinds of files. Run `configure --help' for a list of the directories +you can set and what kinds of files go in them. + + If the package supports it, you can cause programs to be installed +with an extra prefix or suffix on their names by giving `configure' the +option `--program-prefix=PREFIX' or `--program-suffix=SUFFIX'. + +Optional Features +================= + +Some packages pay attention to `--enable-FEATURE' options to +`configure', where FEATURE indicates an optional part of the package. +They may also pay attention to `--with-PACKAGE' options, where PACKAGE +is something like `gnu-as' or `x' (for the X Window System). The +`README' should mention any `--enable-' and `--with-' options that the +package recognizes. + + For packages that use the X Window System, `configure' can usually +find the X include and library files automatically, but if it doesn't, +you can use the `configure' options `--x-includes=DIR' and +`--x-libraries=DIR' to specify their locations. + +Specifying the System Type +========================== + +There may be some features `configure' cannot figure out automatically, +but needs to determine by the type of machine the package will run on. +Usually, assuming the package is built to be run on the _same_ +architectures, `configure' can figure that out, but if it prints a +message saying it cannot guess the machine type, give it the +`--build=TYPE' option. TYPE can either be a short name for the system +type, such as `sun4', or a canonical name which has the form: + + CPU-COMPANY-SYSTEM + +where SYSTEM can have one of these forms: + + OS KERNEL-OS + + See the file `config.sub' for the possible values of each field. If +`config.sub' isn't included in this package, then this package doesn't +need to know the machine type. + + If you are _building_ compiler tools for cross-compiling, you should +use the `--target=TYPE' option to select the type of system they will +produce code for. + + If you want to _use_ a cross compiler, that generates code for a +platform different from the build platform, you should specify the +"host" platform (i.e., that on which the generated programs will +eventually be run) with `--host=TYPE'. + +Sharing Defaults +================ + +If you want to set default values for `configure' scripts to share, you +can create a site shell script called `config.site' that gives default +values for variables like `CC', `cache_file', and `prefix'. +`configure' looks for `PREFIX/share/config.site' if it exists, then +`PREFIX/etc/config.site' if it exists. Or, you can set the +`CONFIG_SITE' environment variable to the location of the site script. +A warning: not all `configure' scripts look for a site script. + +Defining Variables +================== + +Variables not defined in a site shell script can be set in the +environment passed to `configure'. However, some packages may run +configure again during the build, and the customized values of these +variables may be lost. In order to avoid this problem, you should set +them in the `configure' command line, using `VAR=value'. For example: + + ./configure CC=/usr/local2/bin/gcc + +causes the specified `gcc' to be used as the C compiler (unless it is +overridden in the site shell script). Here is a another example: + + /bin/bash ./configure CONFIG_SHELL=/bin/bash + +Here the `CONFIG_SHELL=/bin/bash' operand causes subsequent +configuration-related scripts to be executed by `/bin/bash'. + +`configure' Invocation +====================== + +`configure' recognizes the following options to control how it operates. + +`--help' +`-h' + Print a summary of the options to `configure', and exit. + +`--version' +`-V' + Print the version of Autoconf used to generate the `configure' + script, and exit. + +`--cache-file=FILE' + Enable the cache: use and save the results of the tests in FILE, + traditionally `config.cache'. FILE defaults to `/dev/null' to + disable caching. + +`--config-cache' +`-C' + Alias for `--cache-file=config.cache'. + +`--quiet' +`--silent' +`-q' + Do not print messages saying which checks are being made. To + suppress all normal output, redirect it to `/dev/null' (any error + messages will still be shown). + +`--srcdir=DIR' + Look for the package's source code in directory DIR. Usually + `configure' can determine that directory automatically. + +`configure' also accepts some other, not widely useful, options. Run +`configure --help' for more details. + diff --git a/Makefile.am b/Makefile.am new file mode 100644 index 0000000..b65bd56 --- /dev/null +++ b/Makefile.am @@ -0,0 +1,10 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +################################################ + +SUBDIRS = idl resources src + +salomeinclude_HEADERS = dsccode_config.h + +ACLOCAL_AMFLAGS = -I adm_local + diff --git a/NEWS b/NEWS new file mode 100644 index 0000000..e69de29 diff --git a/README b/README new file mode 100644 index 0000000..ebec0db --- /dev/null +++ b/README @@ -0,0 +1,19 @@ + +Installation +============= +To install this SALOME module use the traditional : + + - ./configure --prefix= --with-kernel= + - make + - make install + +Installation from CVS sources +=============================== +If you get it from CVS, run ./autogen.sh before. + + +Using it +========== +1. Add the module to your SALOME application. +2. Execute the calculation schemas in resources (stream2.xml, calcium4.xml) with + the SALOME module YACS (new supervisor module in SALOME 4.x platform). diff --git a/adm_local/check_Kernel.m4 b/adm_local/check_Kernel.m4 new file mode 100644 index 0000000..414ad00 --- /dev/null +++ b/adm_local/check_Kernel.m4 @@ -0,0 +1,59 @@ +# Check availability of Salome's KERNEL binary distribution +# +# Author : Jerome Roy (CEA, 2003) +# + +AC_DEFUN([AC_CHECK_KERNEL],[ + +AC_CHECKING(for Kernel) + +Kernel_ok=no + +AC_ARG_WITH(kernel, + [ --with-kernel=DIR root directory path of KERNEL build or installation], + KERNEL_DIR="$withval",KERNEL_DIR="") + +if test "x$KERNEL_DIR" = "x" ; then + +# no --with-kernel-dir option used + + if test "x$KERNEL_ROOT_DIR" != "x" ; then + + # KERNEL_ROOT_DIR environment variable defined + KERNEL_DIR=$KERNEL_ROOT_DIR + + else + + # search Kernel binaries in PATH variable + AC_PATH_PROG(TEMP, runSalome) + if test "x$TEMP" != "x" ; then + KERNEL_BIN_DIR=`dirname $TEMP` + KERNEL_DIR=`dirname $KERNEL_BIN_DIR` + fi + + fi +# +fi + +if test -f ${KERNEL_DIR}/bin/salome/runSalome ; then + Kernel_ok=yes + AC_MSG_RESULT(Using Kernel module distribution in ${KERNEL_DIR}) + + if test "x$KERNEL_ROOT_DIR" = "x" ; then + KERNEL_ROOT_DIR=${KERNEL_DIR} + fi + if test "x$KERNEL_SITE_DIR" = "x" ; then + KERNEL_SITE_DIR=${KERNEL_ROOT_DIR} + fi + + AC_SUBST(KERNEL_ROOT_DIR) + AC_SUBST(KERNEL_SITE_DIR) + +else + AC_MSG_WARN("Cannot find compiled Kernel module distribution") +fi + +AC_MSG_RESULT(for Kernel: $Kernel_ok) + +])dnl + diff --git a/adm_local/check_omniorb.m4 b/adm_local/check_omniorb.m4 new file mode 100644 index 0000000..aa7ad58 --- /dev/null +++ b/adm_local/check_omniorb.m4 @@ -0,0 +1,251 @@ + +AC_DEFUN([AC_CHECK_OMNIORB],[ +AC_REQUIRE([AC_PROG_CC])dnl +AC_REQUIRE([AC_PROG_CXX])dnl +AC_REQUIRE([AC_PROG_CPP])dnl +AC_REQUIRE([AC_PROG_CXXCPP])dnl + +AC_CHECKING(for omniORB) +omniORB_ok=yes + +if test "x$PYTHON" = "x" +then + CHECK_PYTHON +fi + +AC_LANG_SAVE +AC_LANG_CPLUSPLUS + +AC_PATH_PROG(OMNIORB_IDL, omniidl) +if test "xOMNIORB_IDL" = "x" +then + omniORB_ok=no + AC_MSG_RESULT(omniORB binaries not in PATH variable) +else + omniORB_ok=yes +fi + +if test "x$omniORB_ok" = "xyes" +then + AC_SUBST(OMNIORB_IDL) + + OMNIORB_BIN=`echo ${OMNIORB_IDL} | sed -e "s,[[^/]]*$,,;s,/$,,;s,^$,.,"` + OMNIORB_ROOT=${OMNIORB_BIN} + # one-level up + OMNIORB_ROOT=`echo ${OMNIORB_ROOT} | sed -e "s,[[^/]]*$,,;s,/$,,;s,^$,.,"` + # + # + if test -d $OMNIORB_ROOT/include ; then + # if $OMNIORB_ROOT/include exists, there are a lot of chance that + # this is omniORB4.x installed via configure && make && make install + OMNIORB_LIB=`echo ${OMNIORB_BIN} | sed -e "s,bin\$,lib,"` + OMNIORB_VERSION=4 + else + # omniORB has been installed old way + OMNIORB_LIB=`echo ${OMNIORB_BIN} | sed -e "s,bin/,lib/,"` + # one-level up again + OMNIORB_ROOT=`echo ${OMNIORB_ROOT} | sed -e "s,[[^/]]*$,,;s,/$,,;s,^$,.,"` + if test -d $OMNIORB_ROOT/include/omniORB4 ; then + OMNIORB_VERSION=4 + else + OMNIORB_VERSION=3 + fi + fi + AC_SUBST(OMNIORB_ROOT) + + OMNIORB_INCLUDES="-I$OMNIORB_ROOT/include -I$OMNIORB_ROOT/include/omniORB${OMNIORB_VERSION} -I$OMNIORB_ROOT/include/COS" + AC_SUBST(OMNIORB_INCLUDES) + + # ENABLE_PTHREADS + + OMNIORB_CXXFLAGS="-DOMNIORB_VERSION=$OMNIORB_VERSION" + case $build_cpu in + sparc*) + OMNIORB_CXXFLAGS="$OMNIORB_CXXFLAGS -D__sparc__" + ;; + *86*) + OMNIORB_CXXFLAGS="$OMNIORB_CXXFLAGS -D__x86__" + ;; + esac + case $build_os in + solaris*) + __OSVERSION__=5 + OMNIORB_CXXFLAGS="$OMNIORB_CXXFLAGS -D__sunos__" + ;; + linux*) + __OSVERSION__=2 + OMNIORB_CXXFLAGS="$OMNIORB_CXXFLAGS -D__linux__" + ;; + esac + AC_SUBST(OMNIORB_CXXFLAGS) + + CPPFLAGS_old=$CPPFLAGS + CPPFLAGS="$CPPFLAGS $OMNIORB_CXXFLAGS $OMNIORB_INCLUDES" + + AC_LANG_CPLUSPLUS + AC_CHECK_HEADER(CORBA.h,omniORB_ok="yes",omniORB_ok="no") + + CPPFLAGS=$CPPFLAGS_old + +fi + +dnl omniORB_ok=yes + +if test "x$omniORB_ok" = "xyes" +then + if test "x$OMNIORB_LIB" = "x/usr/lib" + then + OMNIORB_LDFLAGS="" + else + OMNIORB_LDFLAGS="-L$OMNIORB_LIB" + fi + + LIBS_old=$LIBS + LIBS="$LIBS $OMNIORB_LDFLAGS -lomnithread" + + CXXFLAGS_old=$CXXFLAGS + CXXFLAGS="$CXXFLAGS $OMNIORB_CXXFLAGS $OMNIORB_INCLUDES" + + AC_MSG_CHECKING(whether we can link with omnithreads) + AC_CACHE_VAL(salome_cv_lib_omnithreads,[ + AC_TRY_LINK( +#include +, omni_mutex my_mutex, + eval "salome_cv_lib_omnithreads=yes",eval "salome_cv_lib_omnithreads=no") + ]) + + omniORB_ok="$salome_cv_lib_omnithreads" + if test "x$omniORB_ok" = "xno" + then + AC_MSG_RESULT(omnithreads not found) + else + AC_MSG_RESULT(yes) + fi + + LIBS=$LIBS_old + CXXFLAGS=$CXXFLAGS_old +fi + + +dnl omniORB_ok=yes +if test "x$omniORB_ok" = "xyes" +then + + AC_CHECK_LIB(socket,socket, LIBS="-lsocket $LIBS",,) + AC_CHECK_LIB(nsl,gethostbyname, LIBS="-lnsl $LIBS",,) + + LIBS_old=$LIBS + OMNIORB_LIBS="$OMNIORB_LDFLAGS" + OMNIORB_LIBS="$OMNIORB_LIBS -lomniORB${OMNIORB_VERSION}" + OMNIORB_LIBS="$OMNIORB_LIBS -lomniDynamic${OMNIORB_VERSION}" + OMNIORB_LIBS="$OMNIORB_LIBS -lCOS${OMNIORB_VERSION}" + OMNIORB_LIBS="$OMNIORB_LIBS -lCOSDynamic${OMNIORB_VERSION}" + OMNIORB_LIBS="$OMNIORB_LIBS -lomnithread" + if test $OMNIORB_VERSION = 3 ; then + OMNIORB_LIBS="$OMNIORB_LIBS -ltcpwrapGK" + fi + AC_SUBST(OMNIORB_LIBS) + + LIBS="$OMNIORB_LIBS $LIBS" + CXXFLAGS_old=$CXXFLAGS + CXXFLAGS="$CXXFLAGS $OMNIORB_CXXFLAGS $OMNIORB_INCLUDES" + + AC_MSG_CHECKING(whether we can link with omniORB) + AC_CACHE_VAL(salome_cv_lib_omniorb,[ + AC_TRY_LINK( +#include +, CORBA::ORB_var orb, + eval "salome_cv_lib_omniorb3=yes",eval "salome_cv_lib_omniorb3=no") + ]) + omniORB_ok="$salome_cv_lib_omniorb3" + + omniORB_ok=yes + if test "x$omniORB_ok" = "xno" + then + AC_MSG_RESULT(omniORB library linking failed) + omniORB_ok=no + else + AC_MSG_RESULT(yes) + fi + LIBS="$LIBS_old" + CXXFLAGS=$CXXFLAGS_old +fi + + +if test "x$omniORB_ok" = "xyes" +then + + OMNIORB_IDLCXXFLAGS="-nf -I$OMNIORB_ROOT/idl" + OMNIORB_IDLPYFLAGS="-bpython -I$OMNIORB_ROOT/idl" + AC_SUBST(OMNIORB_IDLCXXFLAGS) + AC_SUBST(OMNIORB_IDLPYFLAGS) + + OMNIORB_IDL_CLN_H=.hh + OMNIORB_IDL_CLN_CXX=SK.cc + OMNIORB_IDL_CLN_OBJ=SK.o + AC_SUBST(OMNIORB_IDL_CLN_H) + AC_SUBST(OMNIORB_IDL_CLN_CXX) + AC_SUBST(OMNIORB_IDL_CLN_OBJ) + IDL_CLN_H=.hh + IDL_CLN_CXX=SK.cc + IDL_CLN_OBJ=SK.o + AC_SUBST(IDL_CLN_H) + AC_SUBST(IDL_CLN_CXX) + AC_SUBST(IDL_CLN_OBJ) + + OMNIORB_IDL_SRV_H=.hh + OMNIORB_IDL_SRV_CXX=SK.cc + OMNIORB_IDL_SRV_OBJ=SK.o + AC_SUBST(OMNIORB_IDL_SRV_H) + AC_SUBST(OMNIORB_IDL_SRV_CXX) + AC_SUBST(OMNIORB_IDL_SRV_OBJ) + IDL_SRV_H=.hh + IDL_SRV_CXX=SK.cc + IDL_SRV_OBJ=SK.o + AC_SUBST(IDL_SRV_H) + AC_SUBST(IDL_SRV_CXX) + AC_SUBST(IDL_SRV_OBJ) + + OMNIORB_IDL_TIE_H= + OMNIORB_IDL_TIE_CXX= + AC_SUBST(OMNIORB_IDL_TIE_H) + AC_SUBST(OMNIORB_IDL_TIE_CXX) + + AC_DEFINE(OMNIORB,,[Presence de omniORB]) + + CORBA_HAVE_POA=1 + AC_DEFINE(CORBA_HAVE_POA,,[POA presence]) + + CORBA_ORB_INIT_HAVE_3_ARGS=1 + AC_DEFINE(CORBA_ORB_INIT_HAVE_3_ARGS,,[?]) + CORBA_ORB_INIT_THIRD_ARG='"omniORB"' + AC_DEFINE(CORBA_ORB_INIT_THIRD_ARG, "omniORB", [?]) + +fi + +omniORBpy_ok=no +if test "x$omniORB_ok" = "xyes" +then + AC_MSG_CHECKING(omniORBpy) + $PYTHON -c "import omniORB" &> /dev/null + if test $? = 0 ; then + AC_MSG_RESULT(yes) + omniORBpy_ok=yes + else + AC_MSG_RESULT(no, check your installation of omniORBpy) + omniORBpy_ok=no + fi +fi + +dnl AC_LANG_RESTORE + +AC_MSG_RESULT(for omniORBpy: $omniORBpy_ok) +AC_MSG_RESULT(for omniORB: $omniORB_ok) + +IDL=${OMNIORB_IDL} +IDLGENFLAGS="-bcxx " +AC_SUBST(IDL) +AC_SUBST(IDLGENFLAGS) + +])dnl +dnl diff --git a/adm_local/make_common_starter.am b/adm_local/make_common_starter.am new file mode 100644 index 0000000..36b1dc3 --- /dev/null +++ b/adm_local/make_common_starter.am @@ -0,0 +1,11 @@ +salomeincludedir = $(prefix)/include/salome +salomelibdir = $(prefix)/lib/salome +salomeidldir = $(prefix)/idl/salome +salomepythondir = $(prefix)/lib/python$(PYTHON_VERSION)/site-packages/salome +salomesharedir= $(prefix)/share/salome/resources + +IDL_INCLUDES = -I$(KERNEL_ROOT_DIR)/idl/salome +KERNEL_LIBS= -L$(KERNEL_ROOT_DIR)/lib/salome -lSalomeContainer -lOpUtil -lSalomeDSCContainer -lSalomeDSCSuperv -lSalomeDatastream -lSalomeDSCSupervBasic -lCalciumC +KERNEL_INCLUDES= -I$(KERNEL_ROOT_DIR)/include/salome $(OMNIORB_INCLUDES) + +C_DEPEND_FLAG=-MM -MG diff --git a/adm_local/python.m4 b/adm_local/python.m4 new file mode 100644 index 0000000..a8013e2 --- /dev/null +++ b/adm_local/python.m4 @@ -0,0 +1,163 @@ +dnl Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +dnl CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +dnl +dnl This library is free software; you can redistribute it and/or +dnl modify it under the terms of the GNU Lesser General Public +dnl License as published by the Free Software Foundation; either +dnl version 2.1 of the License. +dnl +dnl This library is distributed in the hope that it will be useful, +dnl but WITHOUT ANY WARRANTY; without even the implied warranty of +dnl MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +dnl Lesser General Public License for more details. +dnl +dnl You should have received a copy of the GNU Lesser General Public +dnl License along with this library; if not, write to the Free Software +dnl Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +dnl +dnl See http://www.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org +dnl +dnl +dnl +## ------------------------ +## Python file handling +## From Andrew Dalke +## Modified by Marc Tajchman (06/2001) +## ------------------------ + +dnl CHECK_PYTHON([module, classes]) +dnl +dnl Adds support for distributing Python modules or classes. +dnl Python library files distributed as a `module' are installed +dnl under PYTHON_SITE_PACKAGE (eg, ./python1.5/site-package/package-name) +dnl while those distributed as `classes' are installed under PYTHON_SITE +dnl (eg, ./python1.5/site-packages). The default is to install as +dnl a `module'. + +AC_DEFUN([CHECK_PYTHON], + [ + AC_ARG_WITH(python, + [ --with-python=DIR root directory path of python installation ], + [PYTHON="$withval/bin/python" + AC_MSG_RESULT("select python distribution in $withval") + ], [ + AC_PATH_PROG(PYTHON, python) + ]) + + AC_CHECKING([local Python configuration]) + PYTHON_PREFIX=`echo $PYTHON | sed -e "s,[[^/]]*$,,;s,/$,,;s,^$,.,"` + PYTHON_PREFIX=`echo $PYTHON_PREFIX | sed -e "s,[[^/]]*$,,;s,/$,,;s,^$,.,"` + PYTHONHOME=$PYTHON_PREFIX + + AC_SUBST(PYTHON_PREFIX) + AC_SUBST(PYTHONHOME) + + changequote(<<, >>)dnl + PYTHON_VERSION=`$PYTHON -c "import sys; print sys.version[:3]"` + changequote([, ])dnl + AC_SUBST(PYTHON_VERSION) + + PY_MAKEFILE=$PYTHON_PREFIX/lib/python$PYTHON_VERSION/config/Makefile + if test ! -f "$PY_MAKEFILE"; then + AC_MSG_ERROR([*** Couldn't find ${PY_MAKEFILE}. Maybe you are +*** missing the development portion of the python installation]) + fi + + AC_SUBST(PYTHON_INCLUDES) + AC_SUBST(PYTHON_LIBS) + + PYTHON_INCLUDES=-I$PYTHON_PREFIX/include/python$PYTHON_VERSION + PYTHON_LIBS="-L${PYTHON_PREFIX}/lib/python${PYTHON_VERSION}/config -lpython${PYTHON_VERSION}" + PYTHON_LIB=$PYTHON_LIBS + PYTHON_LIBA=$PYTHON_PREFIX/lib/python$PYTHON_VERSION/config/libpython$PYTHON_VERSION.a + + dnl At times (like when building shared libraries) you may want + dnl to know which OS Python thinks this is. + + AC_SUBST(PYTHON_PLATFORM) + PYTHON_PLATFORM=`$PYTHON -c "import sys; print sys.platform"` + + AC_SUBST(PYTHON_SITE) + AC_ARG_WITH(python-site, +[ --with-python-site=DIR Use DIR for installing platform independent + Python site-packages], + +dnl modification : by default, we install python script in salome root tree + +dnl [PYTHON_SITE="$withval" +dnl python_site_given=yes], +dnl [PYTHON_SITE=$PYTHON_PREFIX"/lib/python"$PYTHON_VERSION/site-packages +dnl python_site_given=no]) + +[PYTHON_SITE="$withval" +python_site_given=yes], +[PYTHON_SITE=$prefix"/lib/python"$PYTHON_VERSION/site-packages +python_site_given=no]) + + AC_SUBST(PYTHON_SITE_PACKAGE) + PYTHON_SITE_PACKAGE=$PYTHON_SITE/$PACKAGE + + + dnl Get PYTHON_SITE from --with-python-site-exec or from + dnl --with-python-site or from running Python + + AC_SUBST(PYTHON_SITE_EXEC) + AC_ARG_WITH(python-site-exec, +[ --with-python-site-exec=DIR Use DIR for installing platform dependent + Python site-packages], +[PYTHON_SITE_EXEC="$withval"], +[if test "$python_site_given" = yes; then + PYTHON_SITE_EXEC=$PYTHON_SITE +else + PYTHON_SITE_EXEC=$PYTHON_EXEC_PREFIX"/lib/python"$PYTHON_VERSION/site-packages +fi]) + + dnl Set up the install directory + ifelse($1, classes, +[PYTHON_SITE_INSTALL=$PYTHON_SITE], +[PYTHON_SITE_INSTALL=$PYTHON_SITE_PACKAGE]) + AC_SUBST(PYTHON_SITE_INSTALL) + + dnl Also lets automake think PYTHON means something. + + pythondir=$PYTHON_PREFIX"/lib/python"$PYTHON_VERSION/ + AC_SUBST(pythondir) + + AC_MSG_CHECKING([if we need libdb]) + PY_NEEDOPENDB=`nm $PYTHON_LIBA | grep dbopen | grep U` + if test "x$PY_NEEDOPENDB" != "x"; then + AC_MSG_RESULT(yes) + AC_CHECK_LIB(db,dbopen,PYTHON_LIBS="$PYTHON_LIBS -ldb",db_ok=no) + else + AC_MSG_RESULT(no) + fi + + AC_MSG_CHECKING([if we need libdl]) + PY_NEEDOPENDL=`nm $PYTHON_LIBA | grep dlopen | grep U` + if test "x$PY_NEEDOPENDL" != "x"; then + AC_MSG_RESULT(yes) + AC_CHECK_LIB(dl,dlopen,PYTHON_LIBS="$PYTHON_LIBS -ldl",dl_ok=no) + else + AC_MSG_RESULT(no) + fi + + AC_MSG_CHECKING([if we need libutil]) + PY_NEEDOPENPTY=`nm $PYTHON_LIBA | grep openpty | grep U` + if test "x$PY_NEEDOPENPTY" != "x"; then + AC_MSG_RESULT(yes) + AC_CHECK_LIB(util,openpty,PYTHON_LIBS="$PYTHON_LIBS -lutil",openpty_ok=no) + else + AC_MSG_RESULT(no) + fi + + AC_MSG_CHECKING([if we need tcltk]) + PY_NEEDTCLTK=`nm $PYTHON_LIBA | grep Tcl_Init | grep U` + if test "x$PY_NEEDTCLTK" != "x"; then + AC_MSG_RESULT(yes) + AC_CHECK_LIB(tcl,Tcl_Init,PYTHON_LIBS="$PYTHON_LIBS -ltcl -ltk",tclinit_ok=no) + else + AC_MSG_RESULT(no) + fi + + python_ok=yes + AC_MSG_RESULT(looks good)]) diff --git a/autogen.sh b/autogen.sh new file mode 100755 index 0000000..e4226ae --- /dev/null +++ b/autogen.sh @@ -0,0 +1,12 @@ +#!/bin/sh + +rm -rf autom4te.cache +rm -f aclocal.m4 adm_local/ltmain.sh + +echo "Running aclocal..." ; +aclocal --force -I adm_local || exit 1 +echo "Running autoheader..." ; autoheader --force -I adm_local || exit 1 +echo "Running autoconf..." ; autoconf --force || exit 1 +echo "Running libtoolize..." ; libtoolize --copy --force || exit 1 +echo "Running automake..." ; automake --add-missing --copy --gnu || exit 1 + diff --git a/configure.ac b/configure.ac new file mode 100644 index 0000000..12b9bc9 --- /dev/null +++ b/configure.ac @@ -0,0 +1,37 @@ +AC_INIT(dsccode,0.1) +AC_CONFIG_AUX_DIR(adm_local) +AM_INIT_AUTOMAKE +AC_CONFIG_HEADER(dsccode_config.h) +AC_PROG_LIBTOOL + +AC_PROG_CC +dnl C++ +AC_PROG_CXX + +dnl Check Salome Install +AC_CHECK_KERNEL + +if test "x$Kernel_ok" = "xno"; then + AC_MSG_ERROR([You must define a correct KERNEL_ROOT_DIR !]) +fi + +MACHINE=PCLINUX +AC_SUBST(MACHINE) + +AC_CHECK_OMNIORB + +AC_CONFIG_FILES([ + Makefile + idl/Makefile + resources/Makefile + src/Makefile + src/DSCCODAENG/Makefile + src/DSCCODBENG/Makefile + src/DSCCODCENG/Makefile + src/DSCCODDENG/Makefile + src/NEUTRO/Makefile + src/FLUIDE/Makefile + src/SOLIDE/Makefile + src/INTERPI/Makefile + ]) +AC_OUTPUT diff --git a/idl/DSCCODE.idl b/idl/DSCCODE.idl new file mode 100644 index 0000000..ac8cd13 --- /dev/null +++ b/idl/DSCCODE.idl @@ -0,0 +1,55 @@ +#ifndef _DSCCODE_IDL_ +#define _DSCCODE_IDL_ + +#include "DSC_Engines.idl" +#include "SALOME_Exception.idl" + +module DSCCODE { + + interface DSCCODA: Engines::Superv_Component + { + void prun(in long niter) raises (SALOME::SALOME_Exception); + void trun(in long niter) raises (SALOME::SALOME_Exception); + }; + + interface DSCCODB: Engines::Superv_Component + { + void prun(in long niter) raises (SALOME::SALOME_Exception); + void trun(in long niter) raises (SALOME::SALOME_Exception); + }; + + interface DSCCODC: Engines::Superv_Component + { + void prun(in long niter) raises (SALOME::SALOME_Exception); + void trun(in long niter) raises (SALOME::SALOME_Exception); + }; + + interface DSCCODD: Engines::Superv_Component + { + void prun(in long niter) raises (SALOME::SALOME_Exception); + void trun(in long niter) raises (SALOME::SALOME_Exception); + }; + interface NEUTRO: Engines::Superv_Component + { + void prun() raises (SALOME::SALOME_Exception); + void trun(in double dt) raises (SALOME::SALOME_Exception); + }; + interface FLUIDE: Engines::Superv_Component + { + void prun() raises (SALOME::SALOME_Exception); + void trun(in double dt) raises (SALOME::SALOME_Exception); + }; + interface SOLIDE: Engines::Superv_Component + { + void prun() raises (SALOME::SALOME_Exception); + void trun(in double dt) raises (SALOME::SALOME_Exception); + }; + interface INTERPI: Engines::Superv_Component + { + void prun() raises (SALOME::SALOME_Exception); + void trun() raises (SALOME::SALOME_Exception); + }; +}; + +#endif // _DSCCODE_IDL_ + diff --git a/idl/Makefile.am b/idl/Makefile.am new file mode 100644 index 0000000..d00da26 --- /dev/null +++ b/idl/Makefile.am @@ -0,0 +1,69 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +################################################################# + +BUILT_SOURCES = Calcium_Ports.hh Palm_Ports.hh SALOME_Component.hh SALOME_Ports.hh DSC_Engines.hh SALOME_Exception.hh +IDL_FILES=DSCCODE.idl + +################################################################# + +salomelib_LTLIBRARIES = libDSCCODE.la +salomeidl_DATA = $(IDL_FILES) +salomepython_DATA = DSCCODE_idl.py +libDSCCODE_la_SOURCES = +nodist_libDSCCODE_la_SOURCES = DSCCODESK.cc +libDSCCODE_la_CXXFLAGS = -I. $(KERNEL_INCLUDES) +libDSCCODE_la_LIBADD = $(KERNEL_LIBS) + + +################################################################# +CLEANFILES = *.hh *SK.cc *.py + +clean-local: + rm -rf DSCCODE DSCCODE__POA + +install-data-local: + ${mkinstalldirs} $(DESTDIR)$(salomepythondir) + cp -R DSCCODE DSCCODE__POA $(DESTDIR)$(salomepythondir) + +uninstall-local: + rm -rf $(DESTDIR)$(salomepythondir)/DSCCODE + rm -rf $(DESTDIR)$(salomepythondir)/DSCCODE__POA +################################################################# + +DSC_Engines.hh:$(KERNEL_ROOT_DIR)/idl/salome/DSC_Engines.idl + $(OMNIORB_IDL) $(OMNIORB_IDLCXXFLAGS) -bcxx -I$(KERNEL_ROOT_DIR)/idl/salome $< +SALOME_Ports.hh:$(KERNEL_ROOT_DIR)/idl/salome/SALOME_Ports.idl + $(OMNIORB_IDL) $(OMNIORB_IDLCXXFLAGS) -bcxx -I$(KERNEL_ROOT_DIR)/idl/salome $< +Palm_Ports.hh:$(KERNEL_ROOT_DIR)/idl/salome/Palm_Ports.idl + $(OMNIORB_IDL) $(OMNIORB_IDLCXXFLAGS) -bcxx -I$(KERNEL_ROOT_DIR)/idl/salome $< +Calcium_Ports.hh:$(KERNEL_ROOT_DIR)/idl/salome/Calcium_Ports.idl + $(OMNIORB_IDL) $(OMNIORB_IDLCXXFLAGS) -bcxx -I$(KERNEL_ROOT_DIR)/idl/salome $< + +SALOME_Component.hh:$(KERNEL_ROOT_DIR)/idl/salome/SALOME_Component.idl + $(OMNIORB_IDL) $(OMNIORB_IDLCXXFLAGS) -bcxx -I$(KERNEL_ROOT_DIR)/idl/salome $< +SALOME_Exception.hh:$(KERNEL_ROOT_DIR)/idl/salome/SALOME_Exception.idl + $(OMNIORB_IDL) $(OMNIORB_IDLCXXFLAGS) -bcxx -I$(KERNEL_ROOT_DIR)/idl/salome $< + +%SK.cc %.hh : %.idl + $(IDL) $(IDLGENFLAGS) $(IDL_INCLUDES) $< + +%_idl.py : %.idl + $(IDL) -bpython -I. $(IDL_INCLUDES) $< + +# we use cpp to generate dependencies between idl files. +# option x c tells the preprocessor to consider idl as a c file. +# if an idl is modified, all idl dependencies are rebuilt + +.depidl: $(IDL_FILES) + echo "" > $@ + for dep in $^ dummy; do \ + if [ $$dep != "dummy" ]; then \ + echo Building dependencies for $$dep; \ + $(CPP) $(C_DEPEND_FLAG) -x c -I$(srcdir) $(IDL_INCLUDES) $$dep 2>/dev/null | \ + sed 's/\.o/\SK.cc/' >>$@; \ + fi; \ + done ; + +-include .depidl + diff --git a/resources/DSCCODECatalog.xml b/resources/DSCCODECatalog.xml new file mode 100644 index 0000000..c23c1dc --- /dev/null +++ b/resources/DSCCODECatalog.xml @@ -0,0 +1,330 @@ + + + + + + + + + + + + + + + DSCCODA + DSCCODA + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + DSCCODA.png + 'linux' ~ OS + + + DSCCODA + No comment + + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + DSCCODB + DSCCODB + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + DSCCODB.png + 'linux' ~ OS + + DSCCODB + No comment + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + DSCCODC + DSCCODC + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + DSCCODC.png + 'linux' ~ OS + + DSCCODC + No comment + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + DSCCODD + DSCCODD + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + DSCCODD.png + 'linux' ~ OS + + DSCCODD + No comment + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + SOLIDE + SOLIDE + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + SOLIDE.png + 'linux' ~ OS + + SOLIDE + No comment + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + FLUIDE + FLUIDE + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + FLUIDE.png + 'linux' ~ OS + + FLUIDE + No comment + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + NEUTRO + NEUTRO + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + NEUTRO.png + 'linux' ~ OS + + NEUTRO + No comment + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + + INTERPI + INTERPI + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + INTERPI.png + 'linux' ~ OS + + INTERPI + No comment + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + + INTERPJ + INTERPJ + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + INTERPI.png + 'linux' ~ OS + + INTERPJ + No comment + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + + INTERPK + INTERPK + Data + C. Caremoli + 3.2.0 + EDF - RD + 1 + INTERPI.png + 'linux' ~ OS + + INTERPK + No comment + + + + run + CCar + 1.0 + run + 1 + + + + + + + + + + + + diff --git a/resources/Makefile.am b/resources/Makefile.am new file mode 100644 index 0000000..6922ce4 --- /dev/null +++ b/resources/Makefile.am @@ -0,0 +1,9 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +############################################### +DATA_INST = DSCCODECatalog.xml + +salomeshare_DATA = ${DATA_INST} + +EXTRA_DIST = ${DATA_INST} + diff --git a/resources/calcium4.xml b/resources/calcium4.xml new file mode 100644 index 0000000..cbb9bbe --- /dev/null +++ b/resources/calcium4.xml @@ -0,0 +1,209 @@ + + + + + + + + + + + + + + + + + + + + + + + FLUIDE + prun + + + + + + + SOLIDE + prun + + + + + + + + + NEUTRO + prun + + + + + + + INTERPI + prun + + + + + + crayon tpi + int4 tparoi + + + + int4 tpar + canal tpi + + + + canal tfi + crayon tfi + + + + crayon tempi + comb tempi + + + + comb puissi + crayon puissi + + + + crayon iconv + canal iconv + + + + crayon iconv + comb iconv + + + + + + + a.canal + trun + + + + + + + + a.crayon + trun + + + + + + + + + + a.comb + trun + + + + + + a.int4 + trun + + + + + INTERPJ + trun + + + + + + INTERPK + trun + + + + + + canal rfluide + crayon rext + + + + canal tfluide + crayon text + + + + crayon rparoi + canal rparoi + + + + crayon tparoi + int1 tparoit + + + + int1 tpart + canal tparoi + + + + crayon temperature + int2 tparoit + + + + int2 tpart + comb temperature + + + + comb puissance + int3 tparoit + + + + int3 tpart + crayon puissa + + + + + + a b + + + + + + b.canal dt + 0.8 + + + b.crayon dt + 0.8 + + + b.comb dt + 0.8 + + + + + + diff --git a/resources/stream2.xml b/resources/stream2.xml new file mode 100644 index 0000000..6b1c733 --- /dev/null +++ b/resources/stream2.xml @@ -0,0 +1,51 @@ + + + + + + + + + + + + DSCCODA + prun + + + + + + DSCCODB + prun + + + + + + node0 node1 + node0 node2 + + + node0n + node1 niter + + + node0n + node2 niter + + + + + + node1 prun_out_port + node2 prun_in_port + + + node2 prun_out_port + node1 prun_in_port + + + diff --git a/src/DSCCODAENG/DSCCODAENG.cxx b/src/DSCCODAENG/DSCCODAENG.cxx new file mode 100755 index 0000000..2fe4c92 --- /dev/null +++ b/src/DSCCODAENG/DSCCODAENG.cxx @@ -0,0 +1,141 @@ +using namespace std; +#include "DSCCODAENG.hxx" +#include +#include + +//! Constructor for component "DSCCODA" instance +/*! + * + */ +DSCCODA_i::DSCCODA_i(CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + : Superv_Component_i(orb, poa, contId, instanceName, interfaceName) +{ + cout << "create component" << endl; + _thisObj = this ; + _id = _poa->activate_object(_thisObj); + _in_port=NULL; + _out_port=NULL; +} + +//! Destructor for component "DSCCODA" instance +DSCCODA_i::~DSCCODA_i() +{ +} + +//! Register datastream ports for a component service given its name +/*! + * \param service_name : service name + * \return true if port registering succeeded, false if not + */ +CORBA::Boolean +DSCCODA_i::init_service(const char * service_name) { + CORBA::Boolean rtn = false; + string s_name(service_name); + if (s_name == "prun") + { + try + { + std::cerr << "DSCCODA: prun: uses " << std::endl; + add_port("BASIC_short", "uses", "prun_out_port"); + std::cerr << "DSCCODA: prun: provides " << std::endl; + add_port("BASIC_short", "provides", "prun_in_port"); + } + //catch(const PortAlreadyDefined& ex) + catch(PortAlreadyDefined ex) + { + std::cerr << "DSCCODA: deja defini" << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "DSCCODA: exception inconnue" << std::endl; + } + + rtn = true; + } + if (s_name == "trun") + { + try + { + std::cerr << "DSCCODA: trun: uses " << std::endl; + add_port("BASIC_short", "uses", "trun_out_port"); + std::cerr << "DSCCODA: trun: provides " << std::endl; + add_port("BASIC_short", "provides", "trun_in_port"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "DSCCODA: deja defini" << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "DSCCODA: exception inconnue" << std::endl; + } + rtn = true; + } + return rtn; +} + +void DSCCODA_i::prun(CORBA::Long niter) +{ + cout << "DSCCODA_i::prun" << endl; + uses_port * temp = NULL; + get_port(temp, "prun_out_port"); + _out_port = dynamic_cast(temp); + //get_port((uses_port*&)_out_port, "run_out_port"); + + provides_port * temp2 = NULL; + get_port(temp2, "prun_in_port"); + _in_port = dynamic_cast(temp2); + + int j=0; + _out_port->put(j); + for(int i=0;iget(); + cout << "DSCCODA Got: " << j << endl; + usleep(100000); + j=j+5; + _out_port->put(j); + } +} +void DSCCODA_i::trun(CORBA::Long niter) +{ + cout << "DSCCODA_i::trun" << endl; + uses_port * temp = NULL; + get_port(temp, "trun_out_port"); + _out_port = dynamic_cast(temp); + //get_port((uses_port*&)_out_port, "run_out_port"); + + provides_port * temp2 = NULL; + get_port(temp2, "trun_in_port"); + _in_port = dynamic_cast(temp2); + + int j=0; + for(int i=0;iget(); + cout << "DSCCODA Got: " << j << endl; + usleep(100000); + j=j+5; + _out_port->put(j); + } +} + +extern "C" +{ + PortableServer::ObjectId * DSCCODAEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + MESSAGE("PortableServer::ObjectId * DSCCODAEngine_factory()"); + DSCCODA_i * myEngine = new DSCCODA_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } +} diff --git a/src/DSCCODAENG/DSCCODAENG.hxx b/src/DSCCODAENG/DSCCODAENG.hxx new file mode 100644 index 0000000..e534208 --- /dev/null +++ b/src/DSCCODAENG/DSCCODAENG.hxx @@ -0,0 +1,34 @@ +#ifndef _DSCCODAENG_HXX_ +#define _DSCCODAENG_HXX_ + +#include "Superv_Component_i.hxx" +#include "DSCCODE.hh" +#include "data_short_port_uses.hxx" +#include "data_short_port_provides.hxx" + +class DSCCODA_i: + public virtual POA_DSCCODE::DSCCODA, + public virtual Superv_Component_i +{ + public: + DSCCODA_i(CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, const char *interfaceName); + virtual ~DSCCODA_i(); + CORBA::Boolean init_service(const char * service_name); + void prun(CORBA::Long niter); + void trun(CORBA::Long niter); + private : + data_short_port_uses * _out_port; + data_short_port_provides * _in_port; +}; + +extern "C" +{ + PortableServer::ObjectId * DSCCODAEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName); +} +#endif diff --git a/src/DSCCODAENG/Makefile.am b/src/DSCCODAENG/Makefile.am new file mode 100644 index 0000000..cd29b71 --- /dev/null +++ b/src/DSCCODAENG/Makefile.am @@ -0,0 +1,9 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +salomelib_LTLIBRARIES = libDSCCODAEngine.la +libDSCCODAEngine_la_SOURCES = DSCCODAENG.cxx +nodist_libDSCCODAEngine_la_SOURCES = +libDSCCODAEngine_la_CXXFLAGS = -I$(top_builddir)/idl $(KERNEL_INCLUDES) +libDSCCODAEngine_la_LIBADD = -L$(top_builddir)/idl -lDSCCODE +salomeinclude_HEADERS = DSCCODAENG.hxx + diff --git a/src/DSCCODBENG/DSCCODBENG.cxx b/src/DSCCODBENG/DSCCODBENG.cxx new file mode 100755 index 0000000..35f729a --- /dev/null +++ b/src/DSCCODBENG/DSCCODBENG.cxx @@ -0,0 +1,144 @@ +using namespace std; +#include "DSCCODBENG.hxx" +#include +#include + +//! Constructor for component "DSCCODB" instance +/*! + * + */ +DSCCODB_i::DSCCODB_i(CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + : Superv_Component_i(orb, poa, contId, instanceName, interfaceName) +{ + cout << "create component" << endl; + _thisObj = this ; + _id = _poa->activate_object(_thisObj); + _in_port=NULL; + _out_port=NULL; +} + +//! Destructor for component "DSCCODB" instance +DSCCODB_i::~DSCCODB_i() +{ +} + +//! Register datastream ports for a component service given its name +/*! + * \param service_name : service name + * \return true if port registering succeeded, false if not + */ +CORBA::Boolean +DSCCODB_i::init_service(const char * service_name) { + CORBA::Boolean rtn = false; + string s_name(service_name); + if (s_name == "prun") + { + try + { + std::cerr << "DSCCODB: prun: uses " << std::endl; + add_port("BASIC_short", "uses", "prun_out_port"); + std::cerr << "DSCCODB: prun: provides " << std::endl; + add_port("BASIC_short", "provides", "prun_in_port"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "DSCCODB: prun: deja defini " << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "DSCCODB: exception inconnue" << std::endl; + } + rtn = true; + } + if (s_name == "trun") + { + try + { + std::cerr << "DSCCODB: trun: uses " << std::endl; + add_port("BASIC_short", "uses", "trun_out_port"); + std::cerr << "DSCCODB: trun: provides " << std::endl; + add_port("BASIC_short", "provides", "trun_in_port"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "DSCCODB: trun: deja defini " << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "DSCCODB: exception inconnue" << std::endl; + } + rtn = true; + } + return rtn; +} + +//! Main method for service "prun" +void DSCCODB_i::prun(CORBA::Long niter) +{ + cout << "DSCCODB_i::prun" << endl; + uses_port * temp = NULL; + get_port(temp, "prun_out_port"); + _out_port = dynamic_cast(temp); + //get_port((uses_port*&)_out_port, "prun_out_port"); + + provides_port * temp2 = NULL; + get_port(temp2, "prun_in_port"); + _in_port = dynamic_cast(temp2); + + int j=0; + _out_port->put(j); + for(int i=0;iget(); + cout << "DSCCODB Got: " << j << endl; + usleep(100000); + j=2*j; + _out_port->put(j); + } +} + +//! Main method for service "trun" +void DSCCODB_i::trun(CORBA::Long niter) +{ + cout << "DSCCODB_i::trun" << endl; + uses_port * temp = NULL; + get_port(temp, "trun_out_port"); + _out_port = dynamic_cast(temp); + //get_port((uses_port*&)_out_port, "trun_out_port"); + + provides_port * temp2 = NULL; + get_port(temp2, "trun_in_port"); + _in_port = dynamic_cast(temp2); + + int j=0; + for(int i=0;iget(); + cout << "DSCCODB Got: " << j << endl; + usleep(100000); + j=2*j; + _out_port->put(j); + } +} + + +extern "C" +{ +//! Factory for component "DSCCODB" + PortableServer::ObjectId * DSCCODBEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + MESSAGE("PortableServer::ObjectId * DSCCODBEngine_factory()"); + DSCCODB_i * myEngine = new DSCCODB_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } +} diff --git a/src/DSCCODBENG/DSCCODBENG.hxx b/src/DSCCODBENG/DSCCODBENG.hxx new file mode 100644 index 0000000..0f56e99 --- /dev/null +++ b/src/DSCCODBENG/DSCCODBENG.hxx @@ -0,0 +1,34 @@ +#ifndef _DSCCODBENG_HXX_ +#define _DSCCODBENG_HXX_ + +#include "Superv_Component_i.hxx" +#include "DSCCODE.hh" +#include "data_short_port_uses.hxx" +#include "data_short_port_provides.hxx" + +class DSCCODB_i: + public virtual POA_DSCCODE::DSCCODB, + public virtual Superv_Component_i +{ + public: + DSCCODB_i(CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, const char *interfaceName); + virtual ~DSCCODB_i(); + CORBA::Boolean init_service(const char * service_name); + void prun(CORBA::Long niter); + void trun(CORBA::Long niter); + private: + data_short_port_provides * _in_port; + data_short_port_uses * _out_port; +}; + +extern "C" +{ + PortableServer::ObjectId * DSCCODBEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName); +} +#endif diff --git a/src/DSCCODBENG/Makefile.am b/src/DSCCODBENG/Makefile.am new file mode 100644 index 0000000..e5cbf18 --- /dev/null +++ b/src/DSCCODBENG/Makefile.am @@ -0,0 +1,9 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +salomelib_LTLIBRARIES = libDSCCODBEngine.la +salomeinclude_HEADERS = DSCCODBENG.hxx +libDSCCODBEngine_la_SOURCES = DSCCODBENG.cxx +nodist_libDSCCODBEngine_la_SOURCES = +libDSCCODBEngine_la_CXXFLAGS = -I$(top_builddir)/idl $(KERNEL_INCLUDES) +libDSCCODBEngine_la_LIBADD = -L$(top_builddir)/idl -lDSCCODE + diff --git a/src/DSCCODCENG/DSCCODCENG.cxx b/src/DSCCODCENG/DSCCODCENG.cxx new file mode 100755 index 0000000..7142c53 --- /dev/null +++ b/src/DSCCODCENG/DSCCODCENG.cxx @@ -0,0 +1,162 @@ +using namespace std; +#include "DSCCODCENG.hxx" +#include +#include +#include "CalciumInterface.hxx" +#include + +//! Constructor for component "DSCCODC" instance +/*! + * + */ +DSCCODC_i::DSCCODC_i(CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + : Superv_Component_i(orb, poa, contId, instanceName, interfaceName) +{ + cerr << "create component" << endl; + _thisObj = this ; + _id = _poa->activate_object(_thisObj); +} + +//! Destructor for component "DSCCODC" instance +DSCCODC_i::~DSCCODC_i() +{ +} + +//! Register datastream ports for a component service given its name +/*! + * \param service_name : service name + * \return true if port registering succeeded, false if not + */ +CORBA::Boolean +DSCCODC_i::init_service(const char * service_name) { + CORBA::Boolean rtn = false; + string s_name(service_name); + if (s_name == "prun") + { + try + { + std::cerr << "DSCCODC: prun: " << std::endl; + + //initialization CALCIUM ports + calcium_integer_port_provides * p1; + /* ETP_EN T IN ENTIER */ + p1 = add_port("CALCIUM_integer","provides","ETP_EN"); + p1->setDependencyType(CalciumTypes::TIME_DEPENDENCY); + /* STP_EN T OUT ENTIER */ + add_port("CALCIUM_integer","uses","STP_EN"); + //end of initialization CALCIUM ports + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "DSCCODC: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "DSCCODC: unknown exception" << std::endl; + } + + rtn = true; + } + if (s_name == "trun") + { + try + { + std::cerr << "DSCCODC: trun: " << std::endl; + calcium_real_port_provides * p1; + /* ETP_RE T IN REEL */ + p1 = add_port("CALCIUM_real","provides","ETP_RE"); + p1->setDependencyType(CalciumTypes::TIME_DEPENDENCY); + /* STP_RE T OUT REEL */ + add_port("CALCIUM_real","uses","STP_RE"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "DSCCODC: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "DSCCODC: unknown exception" << std::endl; + } + rtn = true; + } + return rtn; +} + +void DSCCODC_i::prun(CORBA::Long niter) +{ + cerr << "DSCCODC_i::prun" << endl; + Superv_Component_i * component = dynamic_cast(this); + std::cerr << "Valeur de component : " << component << std::endl; + std::cerr << "Valeur de this : " << this << std::endl; + char nom_instance[INSTANCE_LEN]; + int SDATA_EN [3]; // buffer + int EDATA_EN [3]; // buffer + + int info = cp_cd(component,nom_instance); + + float ti_re = 0.0; // step time + int i = 0; // not used + SDATA_EN[0] = 3; + SDATA_EN[1] = 1; + SDATA_EN[2] = 2; + info = cp_een(component,CP_TEMPS,ti_re,i,"STP_EN",2,SDATA_EN); + + float tf_re = 1.0; // step start time + int max=3; // max size expected + int n; // real size received + + info = cp_len(component,CP_TEMPS,&ti_re,&tf_re,&i,"ETP_EN",max,&n,EDATA_EN); + + for (int i = 0; i < n; ++i) + cerr << "seqLongData[" << i << "] = " << EDATA_EN[i] << endl; + + cerr << "end of DSCCODC_i::prun" << endl; +} + +void DSCCODC_i::trun(CORBA::Long niter) +{ + cerr << "DSCCODC_i::trun" << endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + + float SDATA_RE [3], EDATA_RE [3]; + SDATA_RE[0] = 3.; + SDATA_RE[1] = 1.; + SDATA_RE[2] = 2.; + + float ti_re = 0.0; // time + int i = 0; // not used + + info = cp_ere(component,CP_TEMPS,ti_re,i,"STP_RE",2,SDATA_RE); + cerr << "apres cp_ere: " << info << endl; + + int max=3; // max size expected + int n=2; // real size received + float tf_re = 1.0; // time + info = cp_lre(component,CP_TEMPS,&ti_re,&tf_re,&i,"ETP_RE",max,&n,EDATA_RE); + + for (int i = 0; i < n; ++i) + cerr << "seqRealData[" << i << "] = " << EDATA_RE[i] << endl; + cerr << "end of DSCCODC_i::trun" << endl; +} + +extern "C" +{ + PortableServer::ObjectId * DSCCODCEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + MESSAGE("PortableServer::ObjectId * DSCCODCEngine_factory()"); + DSCCODC_i * myEngine = new DSCCODC_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } +} diff --git a/src/DSCCODCENG/DSCCODCENG.hxx b/src/DSCCODCENG/DSCCODCENG.hxx new file mode 100644 index 0000000..23a9125 --- /dev/null +++ b/src/DSCCODCENG/DSCCODCENG.hxx @@ -0,0 +1,29 @@ +#ifndef _DSCCODCENG_HXX_ +#define _DSCCODCENG_HXX_ + +#include "Superv_Component_i.hxx" +#include "DSCCODE.hh" + +class DSCCODC_i: + public virtual POA_DSCCODE::DSCCODC, + public virtual Superv_Component_i +{ + public: + DSCCODC_i(CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, const char *interfaceName); + virtual ~DSCCODC_i(); + CORBA::Boolean init_service(const char * service_name); + void prun(CORBA::Long niter); + void trun(CORBA::Long niter); +}; + +extern "C" +{ + PortableServer::ObjectId * DSCCODCEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName); +} +#endif diff --git a/src/DSCCODCENG/Makefile.am b/src/DSCCODCENG/Makefile.am new file mode 100644 index 0000000..71637b5 --- /dev/null +++ b/src/DSCCODCENG/Makefile.am @@ -0,0 +1,9 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +salomelib_LTLIBRARIES = libDSCCODCEngine.la +salomeinclude_HEADERS = DSCCODCENG.hxx +libDSCCODCEngine_la_SOURCES = DSCCODCENG.cxx +nodist_libDSCCODCEngine_la_SOURCES = +libDSCCODCEngine_la_CXXFLAGS = -I$(top_builddir)/idl $(KERNEL_INCLUDES) +libDSCCODCEngine_la_LIBADD = -L$(top_builddir)/idl -lDSCCODE + diff --git a/src/DSCCODDENG/DSCCODDENG.cxx b/src/DSCCODDENG/DSCCODDENG.cxx new file mode 100755 index 0000000..0bbe95c --- /dev/null +++ b/src/DSCCODDENG/DSCCODDENG.cxx @@ -0,0 +1,160 @@ +using namespace std; +#include "DSCCODDENG.hxx" +#include +#include +#include "CalciumInterface.hxx" +#include + +//! Constructor for component "DSCCODD" instance +/*! + * + */ +DSCCODD_i::DSCCODD_i(CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + : Superv_Component_i(orb, poa, contId, instanceName, interfaceName) +{ + cerr << "create component" << endl; + _thisObj = this ; + _id = _poa->activate_object(_thisObj); +} + +//! Destructor for component "DSCCODD" instance +DSCCODD_i::~DSCCODD_i() +{ +} + +//! Register datastream ports for a component service given its name +/*! + * \param service_name : service name + * \return true if port registering succeeded, false if not + */ +CORBA::Boolean +DSCCODD_i::init_service(const char * service_name) { + CORBA::Boolean rtn = false; + string s_name(service_name); + if (s_name == "prun") + { + try + { + std::cerr << "DSCCODD: prun: " << std::endl; + //initialization CALCIUM ports + calcium_integer_port_provides * p1; + /* ETP_EN T IN ENTIER */ + p1 = add_port("CALCIUM_integer","provides","ETP_EN"); + p1->setDependencyType(CalciumTypes::TIME_DEPENDENCY); + /* STP_EN T OUT ENTIER */ + add_port("CALCIUM_integer","uses","STP_EN"); + //end of initialization CALCIUM ports + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "DSCCODD: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "DSCCODD: unknown exception" << std::endl; + } + + rtn = true; + } + if (s_name == "trun") + { + try + { + std::cerr << "DSCCODD: trun: " << std::endl; + calcium_real_port_provides * p1; + /* ETP_RE T IN REEL */ + p1 = add_port("CALCIUM_real","provides","ETP_RE"); + p1->setDependencyType(CalciumTypes::TIME_DEPENDENCY); + /* STP_RE T OUT REEL */ + add_port("CALCIUM_real","uses","STP_RE"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "DSCCODD: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "DSCCODD: unknown exception" << std::endl; + } + rtn = true; + } + return rtn; +} + +void DSCCODD_i::prun(CORBA::Long niter) +{ + cerr << "DSCCODD_i::prun" << endl; + Superv_Component_i * component = dynamic_cast(this); + std::cerr << "Valeur de component : " << component << std::endl; + std::cerr << "Valeur de this : " << this << std::endl; + char nom_instance[INSTANCE_LEN]; + int SDATA_EN [3]; // buffer + int EDATA_EN [3]; // buffer + + int info = cp_cd(component,nom_instance); + + int i = 0; // not used + float ti_re = 0.0; // step time + float tf_re = 1.0; // step start time + int max=3; // max size expected + int n; // real size received + + info = cp_len(component,CP_TEMPS,&ti_re,&tf_re,&i,"ETP_EN",max,&n,EDATA_EN); + + for (int i = 0; i < n; ++i) + { + SDATA_EN[i]=EDATA_EN[i]*2; + cerr << "seqLongData[" << i << "] = " << EDATA_EN[i] << endl; + } + + info = cp_een(component,CP_TEMPS,ti_re,i,"STP_EN",n,SDATA_EN); + cerr << "end of DSCCODD_i::prun" << endl; +} + +void DSCCODD_i::trun(CORBA::Long niter) +{ + cerr << "DSCCODD_i::trun" << endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + + float SDATA_RE [3], EDATA_RE [3]; + + int i = 0; // not used + float ti_re = 0.0; // time + int max=3; // max size expected + int n=2; // real size received + float tf_re = 1.0; // time + + info = cp_lre(component,CP_TEMPS,&ti_re,&tf_re,&i,"ETP_RE",max,&n,EDATA_RE); + + cerr << "apres cp_lre: " << n << endl; + for (int i = 0; i < n; ++i) + { + SDATA_RE[i]=EDATA_RE[i]*2; + cerr << "seqRealData[" << i << "] = " << EDATA_RE[i] << endl; + } + + info = cp_ere(component,CP_TEMPS,ti_re,i,"STP_RE",n,SDATA_RE); + cerr << "end of DSCCODD_i::trun" << endl; +} + +extern "C" +{ + PortableServer::ObjectId * DSCCODDEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + MESSAGE("PortableServer::ObjectId * DSCCODDEngine_factory()"); + DSCCODD_i * myEngine = new DSCCODD_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } +} diff --git a/src/DSCCODDENG/DSCCODDENG.hxx b/src/DSCCODDENG/DSCCODDENG.hxx new file mode 100644 index 0000000..221cc73 --- /dev/null +++ b/src/DSCCODDENG/DSCCODDENG.hxx @@ -0,0 +1,31 @@ +#ifndef _DSCCODDENG_HXX_ +#define _DSCCODDENG_HXX_ + +#include "Superv_Component_i.hxx" +#include "DSCCODE.hh" +#include "data_short_port_uses.hxx" +#include "data_short_port_provides.hxx" + +class DSCCODD_i: + public virtual POA_DSCCODE::DSCCODD, + public virtual Superv_Component_i +{ + public: + DSCCODD_i(CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, const char *interfaceName); + virtual ~DSCCODD_i(); + CORBA::Boolean init_service(const char * service_name); + void prun(CORBA::Long niter); + void trun(CORBA::Long niter); +}; + +extern "C" +{ + PortableServer::ObjectId * DSCCODDEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName); +} +#endif diff --git a/src/DSCCODDENG/Makefile.am b/src/DSCCODDENG/Makefile.am new file mode 100644 index 0000000..6dc5b6f --- /dev/null +++ b/src/DSCCODDENG/Makefile.am @@ -0,0 +1,9 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +salomelib_LTLIBRARIES = libDSCCODDEngine.la +libDSCCODDEngine_la_SOURCES = DSCCODDENG.cxx +nodist_libDSCCODDEngine_la_SOURCES = +libDSCCODDEngine_la_CXXFLAGS = -I$(top_builddir)/idl $(KERNEL_INCLUDES) +libDSCCODDEngine_la_LIBADD = -L$(top_builddir)/idl -lDSCCODE +salomeinclude_HEADERS = DSCCODDENG.hxx + diff --git a/src/FLUIDE/FLUIDE.cxx b/src/FLUIDE/FLUIDE.cxx new file mode 100755 index 0000000..ccd9271 --- /dev/null +++ b/src/FLUIDE/FLUIDE.cxx @@ -0,0 +1,154 @@ +#include "FLUIDE.hxx" +#include +#include + +#include +#include + +using namespace std; + +extern "C" void transit_(void *compo,double *dt); +extern "C" void perma_(void *compo); + +//! Constructor for component "FLUIDE" instance +/*! + * + */ +FLUIDE_i::FLUIDE_i(CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + : Superv_Component_i(orb, poa, contId, instanceName, interfaceName) +{ + cerr << "create component" << endl; + _thisObj = this ; + _id = _poa->activate_object(_thisObj); +} + +//! Destructor for component "FLUIDE" instance +FLUIDE_i::~FLUIDE_i() +{ +} + +//! Register datastream ports for a component service given its name +/*! + * \param service_name : service name + * \return true if port registering succeeded, false if not + */ +CORBA::Boolean +FLUIDE_i::init_service(const char * service_name) { + CORBA::Boolean rtn = false; + string s_name(service_name); + if (s_name == "prun") + { + try + { + std::cerr << "FLUIDE: prun: " << std::endl; + //initialization CALCIUM ports IN + create_calcium_port(this,"tpi","CALCIUM_real","IN","I"); + create_calcium_port(this,"iconv","CALCIUM_integer","IN","I"); + //initialization CALCIUM ports OUT + create_calcium_port(this,"tfi","CALCIUM_real","OUT","I"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "FLUIDE: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "FLUIDE: unknown exception" << std::endl; + } + + rtn = true; + } + if (s_name == "trun") + { + try + { + std::cerr << "FLUIDE: trun: " << std::endl; + //initialization CALCIUM ports IN + create_calcium_port(this,"tparoi","CALCIUM_real","IN","T"); + create_calcium_port(this,"rparoi","CALCIUM_real","IN","T"); + //initialization CALCIUM ports OUT + create_calcium_port(this,"tfluide","CALCIUM_real","OUT","T"); + create_calcium_port(this,"rfluide","CALCIUM_real","OUT","T"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "FLUIDE: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "FLUIDE: unknown exception" << std::endl; + } + rtn = true; + } + return rtn; +} + +void FLUIDE_i::prun() +{ + std::cerr << "FLUIDE_i::prun" << std::endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + try + { + perma_(&component); + //to do or not to do ??? + cp_fin(component,CP_ARRET); + } + catch ( const CalciumException & ex) + { + std::cerr << ex.what() << std::endl; + cp_fin(component,CP_ARRET); + } + catch (...) + { + std::cerr << "unexpected exception" << std::endl; + cp_fin(component,CP_ARRET); + } + std::cerr << "end of FLUIDE_i::prun" << std::endl; +} + +void FLUIDE_i::trun(CORBA::Double dt) +{ + std::cerr << "FLUIDE_i::trun" << std::endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + try + { + transit_(&component,&dt); + //to do or not to do ??? + cp_fin(component,CP_ARRET); + } + catch ( const CalciumException & ex) + { + std::cerr << ex.what() << std::endl; + cp_fin(component,CP_ARRET); + } + catch (...) + { + std::cerr << "unexpected exception" << std::endl; + cp_fin(component,CP_ARRET); + } + std::cerr << "end of FLUIDE_i::trun" << std::endl; +} + +extern "C" +{ + PortableServer::ObjectId * FLUIDEEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + MESSAGE("PortableServer::ObjectId * FLUIDEEngine_factory()"); + FLUIDE_i * myEngine = new FLUIDE_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } +} diff --git a/src/FLUIDE/FLUIDE.hxx b/src/FLUIDE/FLUIDE.hxx new file mode 100644 index 0000000..e8740d2 --- /dev/null +++ b/src/FLUIDE/FLUIDE.hxx @@ -0,0 +1,29 @@ +#ifndef _FLUIDE_HXX_ +#define _FLUIDE_HXX_ + +#include "Superv_Component_i.hxx" +#include "DSCCODE.hh" + +class FLUIDE_i: + public virtual POA_DSCCODE::FLUIDE, + public virtual Superv_Component_i +{ + public: + FLUIDE_i(CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, const char *interfaceName); + virtual ~FLUIDE_i(); + CORBA::Boolean init_service(const char * service_name); + void prun(); + void trun(CORBA::Double dt); +}; + +extern "C" +{ + PortableServer::ObjectId * FLUIDEEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName); +} +#endif diff --git a/src/FLUIDE/Makefile.am b/src/FLUIDE/Makefile.am new file mode 100644 index 0000000..f3a152e --- /dev/null +++ b/src/FLUIDE/Makefile.am @@ -0,0 +1,12 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +AM_CFLAGS=$(KERNEL_INCLUDES) -fexceptions + +salomelib_LTLIBRARIES = libFLUIDEEngine.la +libFLUIDEEngine_la_SOURCES = FLUIDE.cxx fluide.f +nodist_libFLUIDEEngine_la_SOURCES = +libFLUIDEEngine_la_CXXFLAGS = -I$(top_builddir)/idl $(KERNEL_INCLUDES) +libFLUIDEEngine_la_FFLAGS = $(KERNEL_INCLUDES) -fexceptions +libFLUIDEEngine_la_LIBADD = -L$(top_builddir)/idl -lDSCCODE -lg2c +salomeinclude_HEADERS = FLUIDE.hxx + diff --git a/src/FLUIDE/fluide.f b/src/FLUIDE/fluide.f new file mode 100644 index 0000000..20c3efb --- /dev/null +++ b/src/FLUIDE/fluide.f @@ -0,0 +1,161 @@ + subroutine transit(compo,dt0) + + include "calcium.hf" + dimension tn(7),t(7),rf(7) + dimension tparoi(7),rparoi(7) + dimension maille(2),tflu(2,7) + double precision dt0 + integer compo + + dt=4.5*dt0 + + r=1. + ro=1. + rext=0.5 + te=10. + deb=10. + + maille(1)=2 + maille(2)=7 + + ti=0. + + do i=1,7 + t(i)=100. + tn(i)=100. + rf(i)=0.5 + enddo +c +c initialisation de temperature fluide a t=0 +c + CALL cpeRE(compo,CP_TEMPS, ti, npas+1, 'tfluide', 6, + & t(2) , info) + IF( info.NE. CPOK )GO TO 9000 +c +c initialisation de la resistance thermique fluide a t=0 +c + CALL cpeRE(compo,CP_TEMPS, ti, npas+1, 'rfluide', 6, + & rf(2) , info) + IF( info.NE. CPOK )GO TO 9000 +c +c boucle temporelle jusqu'a 100. +c + do while( ti.lT.100. ) + tf=ti+dt +c +c lecture de la temperature de paroi entre ti et tf +c + CALL cplRE(compo,CP_TEMPS,ti, tf, npas+1,'tparoi', 6, + & nval, tparoi(2), info) + IF( info.NE. CPOK )GO TO 9000 +c +c lecture de la resistance solide de bord entre ti et tf +c + CALL cplRE(compo,CP_TEMPS,ti, tf, npas+1, 'rparoi', 6, + & nval, rparoi(2), info) + IF( info.NE. CPOK )GO TO 9000 +c +c calcul de la temperature a tf +c + do i=2,7 + smb=ro/dt*tn(i)+deb*t(i-1)+tparoi(i)/(rparoi(i)+rf(i)) + t(i)=smb/(ro/dt+deb+1./(rparoi(i)+rf(i))) + enddo + + write(6,*)'FLUID:temps=',tf,' temperature de sortie canal=',t(7) + call flush(6) +c +c ecriture de la temperature fluide a tf +c + CALL cpeRE(compo,CP_TEMPS, tf, npas+1, 'tfluide', 6, + & t(2) , info) + IF( info.NE.CPOK )GO TO 9000 +c +c ecriture de la resistance thermique fluide a tf +c + CALL cpeRE(compo,CP_TEMPS, tf, npas+1, 'rfluide', 6, + & rf(2) , info) + IF( info.NE.CPOK )GO TO 9000 + do i=1,7 + tflu(1,i)=t(i) + tflu(2,i)=t(i) + enddo + + do i=1,7 + tn(i)=t(i) + enddo + + ti=tf + + enddo + +9000 continue + CALL cpfin(compo,CP_ARRET,info) + end + + subroutine perma(compo) + include "calcium.hf" + dimension t(7),rf(7) + dimension tparoi(7),rparoi(7) + integer compo + + deb=10. + + do i=1,7 + t(i)=100. + rf(i)=0.5 + rparoi(i)=0.5 + enddo +c +c initialisation de temperature fluide a i=0 +c + iter=0 + iconv=0 + CALL cpeRE(compo,CP_ITERATION, ti, iter , 'tfi', 6, + & t(2) , info) + IF( info.NE. CPOK )GO TO 9000 +c +c boucle temporelle jusqu'a 100. +c + do while( iconv .EQ. 0) +c +c lecture de la temperature de paroi iteration iter +c + CALL cplRE(compo,CP_ITERATION,ti, tf, iter , 'tpi', 6, + & nval, tparoi(2), info) + IF( info.NE. CPOK )GO TO 9000 +c +c calcul de la temperature +c + do i=2,7 + smb=deb*t(i-1)+tparoi(i)/(rparoi(i)+rf(i)) + t(i)=smb/(deb+1./(rparoi(i)+rf(i))) + enddo +c +c ecriture de la temperature fluide a iter+1 +c + CALL cpeRE(compo,CP_ITERATION,ti,iter+1, 'tfi', 6, + & t(2) , info) + IF( info.NE. CPOK )GO TO 9000 + + iter=iter+1 + write(6,*)'iter = ',iter,' temperature de sortie canal = ',t(7) +c +c lecture du flag de convergence iconv +c + CALL cplEN(compo,CP_ITERATION,ti, tf, iter , 'iconv', 1, + & nval, iconv , info) + write(6,*)"info:",info + write(6,*)"FLUIDE:",iter,iconv + call flush(6) + + IF( info.NE. CPOK )GO TO 9000 + + if(iconv.eq.1)go to 9000 + + enddo + +9000 continue + CALL cpfin(compo,CP_ARRET,info) + end + diff --git a/src/INTERPI/INTERPI.cxx b/src/INTERPI/INTERPI.cxx new file mode 100755 index 0000000..9fbe4a9 --- /dev/null +++ b/src/INTERPI/INTERPI.cxx @@ -0,0 +1,175 @@ +using namespace std; +#include "INTERPI.hxx" +#include +#include + +#include +#include + +extern "C" void transit_(void *compo); +extern "C" void perma_(void *compo); + +//! Constructor for component "INTERPI" instance +/*! + * + */ +INTERPI_i::INTERPI_i(CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + : Superv_Component_i(orb, poa, contId, instanceName, interfaceName) +{ + cerr << "create component" << endl; + _thisObj = this ; + _id = _poa->activate_object(_thisObj); +} + +//! Destructor for component "INTERPI" instance +INTERPI_i::~INTERPI_i() +{ +} + +//! Register datastream ports for a component service given its name +/*! + * \param service_name : service name + * \return true if port registering succeeded, false if not + */ +CORBA::Boolean +INTERPI_i::init_service(const char * service_name) { + CORBA::Boolean rtn = false; + string s_name(service_name); + if (s_name == "prun") + { + try + { + std::cerr << "INTERPI: prun: " << std::endl; + //initialization CALCIUM ports IN + calcium_real_port_provides * p1; + p1 = add_port("CALCIUM_real","provides","tparoi"); + p1->setDependencyType(CalciumTypes::ITERATION_DEPENDENCY); + //initialization CALCIUM ports OUT + add_port("CALCIUM_real","uses","tpar"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "INTERPI: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "INTERPI: unknown exception" << std::endl; + } + + rtn = true; + } + if (s_name == "trun") + { + try + { + std::cerr << "INTERPI: trun: " << std::endl; + //initialization CALCIUM ports IN + calcium_real_port_provides * p1; + p1 = add_port("CALCIUM_real","provides","tparoit"); + p1->setDependencyType(CalciumTypes::TIME_DEPENDENCY); + //initialization CALCIUM ports OUT + add_port("CALCIUM_real","uses","tpart"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "INTERPI: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "INTERPI: unknown exception" << std::endl; + } + rtn = true; + } + return rtn; +} + +void INTERPI_i::prun() +{ + std::cerr << "INTERPI_i::prun" << std::endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + try + { + perma_(&component); + //to do or not to do ??? + cp_fin(component,CP_ARRET); + } + catch ( const CalciumException & ex) + { + std::cerr << ex.what() << std::endl; + cp_fin(component,CP_ARRET); + } + catch (...) + { + std::cerr << "unexpected exception" << std::endl; + cp_fin(component,CP_ARRET); + } + std::cerr << "end of INTERPI_i::prun" << std::endl; +} + +void INTERPI_i::trun() +{ + std::cerr << "INTERPI_i::trun" << std::endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + try + { + transit_(&component); + //to do or not to do ??? + cp_fin(component,CP_ARRET); + } + catch ( const CalciumException & ex) + { + std::cerr << ex.what() << std::endl; + cp_fin(component,CP_ARRET); + } + catch (...) + { + std::cerr << "unexpected exception" << std::endl; + cp_fin(component,CP_ARRET); + } + std::cerr << "end of INTERPI_i::trun" << std::endl; +} + +extern "C" +{ + PortableServer::ObjectId * INTERPIEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + std::cerr << "INTERPIEngine_factory" << std::endl; + INTERPI_i * myEngine = new INTERPI_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } + + PortableServer::ObjectId * INTERPJEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + std::cerr << "INTERPJEngine_factory" << std::endl; + INTERPI_i * myEngine = new INTERPI_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } + PortableServer::ObjectId * INTERPKEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + std::cerr << "INTERPKEngine_factory" << std::endl; + INTERPI_i * myEngine = new INTERPI_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } +} diff --git a/src/INTERPI/INTERPI.hxx b/src/INTERPI/INTERPI.hxx new file mode 100644 index 0000000..760789d --- /dev/null +++ b/src/INTERPI/INTERPI.hxx @@ -0,0 +1,27 @@ +#ifndef _INTERPI_HXX_ +#define _INTERPI_HXX_ + +#include "Superv_Component_i.hxx" +#include "DSCCODE.hh" + +class INTERPI_i: + public virtual POA_DSCCODE::INTERPI, + public virtual Superv_Component_i +{ + public: + INTERPI_i(CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, const char *interfaceName); + virtual ~INTERPI_i(); + CORBA::Boolean init_service(const char * service_name); + void prun(); + void trun(); +}; + +extern "C" +{ + PortableServer::ObjectId * INTERPIEngine_factory( CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, PortableServer::ObjectId * contId, const char *instanceName, const char *interfaceName); + PortableServer::ObjectId * INTERPJEngine_factory( CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, PortableServer::ObjectId * contId, const char *instanceName, const char *interfaceName); + PortableServer::ObjectId * INTERPKEngine_factory( CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, PortableServer::ObjectId * contId, const char *instanceName, const char *interfaceName); +} +#endif diff --git a/src/INTERPI/Makefile.am b/src/INTERPI/Makefile.am new file mode 100644 index 0000000..1cf79db --- /dev/null +++ b/src/INTERPI/Makefile.am @@ -0,0 +1,15 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +AM_CFLAGS=$(KERNEL_INCLUDES) -fexceptions +AM_FFLAGS=$(KERNEL_INCLUDES) -fexceptions +AM_CXXFLAGS = -I$(top_builddir)/idl $(KERNEL_INCLUDES) +AM_LDFLAGS = -L$(top_builddir)/idl -lDSCCODE -lg2c + +salomelib_LTLIBRARIES = libINTERPIEngine.la libINTERPJEngine.la libINTERPKEngine.la + +libINTERPIEngine_la_SOURCES = INTERPI.cxx interi.f inter.f +libINTERPJEngine_la_SOURCES = INTERPI.cxx interi.f inter.f +libINTERPKEngine_la_SOURCES = INTERPI.cxx interi.f inter.f + +salomeinclude_HEADERS = INTERPI.hxx + diff --git a/src/INTERPI/inter.f b/src/INTERPI/inter.f new file mode 100644 index 0000000..69479a2 --- /dev/null +++ b/src/INTERPI/inter.f @@ -0,0 +1,45 @@ + subroutine transit(compo) + + include "calcium.hf" + + dimension tparoi(100),tp(100) + character*(INSTANCE_LEN) nom_instance + integer info + integer compo +c +c boucle temporelle infinie +c + do while( .TRUE. ) +c +c lecture de la temperature de paroi a t +c + t=0. + CALL cplRE(compo,CP_SEQUENTIEL,t, t, ii, 'tparoit', 100, + & nval, tparoi, info) + + IF( info.NE.CPOK ) GO TO 9000 +c +c Ici on realise l'interpolation +c + do i=1,nval + tp(i)=tparoi(i) + enddo +c +c ecriture de la temperature de paroi interpolee a t +c + write(6,*)'INTERP:temps=',t + call flush(6) + + CALL cpere(compo,CP_TEMPS, t, ii, 'tpart', nval, + & tp , info) + + IF( info.NE.CPOK ) GO TO 9000 +c if(t .gt. 100.)GO TO 9000 + + enddo + +9000 continue + + CALL cpfin(compo,CP_ARRET, info) + + end diff --git a/src/INTERPI/interi.f b/src/INTERPI/interi.f new file mode 100644 index 0000000..380557e --- /dev/null +++ b/src/INTERPI/interi.f @@ -0,0 +1,39 @@ + subroutine perma(compo) + include "calcium.hf" + + dimension tparoi(100),tp(100) + character*(INSTANCE_LEN) nom_instance + integer info + integer compo + +c +c boucle temporelle infinie +c + do while( .TRUE. ) +c +c lecture de la temperature de paroi a t +c + CALL cplRE(compo,CP_SEQUENTIEL,t, t, ii, 'tparoi', 100, + & nval, tparoi, info) + + IF( info.NE. CPOK ) GO TO 9000 +c +c Ici on realise l'interpolation +c + do i=1,nval + tp(i)=tparoi(i) + enddo +c +c ecriture de la temperature de paroi interpolee a t +c + CALL cpeRE(compo,CP_ITERATION, t, ii, 'tpar', nval, + & tp , info) + + IF( info.NE.CPOK ) GO TO 9000 +c if(ii .GE. 24)GO TO 9000 + + enddo + +9000 continue + CALL cpfin(compo,CP_ARRET, info) + end diff --git a/src/Makefile.am b/src/Makefile.am new file mode 100644 index 0000000..210162b --- /dev/null +++ b/src/Makefile.am @@ -0,0 +1 @@ +SUBDIRS =DSCCODAENG DSCCODBENG DSCCODCENG DSCCODDENG NEUTRO FLUIDE SOLIDE INTERPI diff --git a/src/NEUTRO/Makefile.am b/src/NEUTRO/Makefile.am new file mode 100644 index 0000000..62cb0eb --- /dev/null +++ b/src/NEUTRO/Makefile.am @@ -0,0 +1,12 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +AM_CFLAGS=$(KERNEL_INCLUDES) -fexceptions + +salomelib_LTLIBRARIES = libNEUTROEngine.la +libNEUTROEngine_la_SOURCES = NEUTRO.cxx neutro.f +nodist_libNEUTROEngine_la_SOURCES = +libNEUTROEngine_la_CXXFLAGS = -I$(top_builddir)/idl $(KERNEL_INCLUDES) +libNEUTROEngine_la_FFLAGS = $(KERNEL_INCLUDES) -fexceptions +libNEUTROEngine_la_LIBADD = -L$(top_builddir)/idl -lDSCCODE -lg2c +salomeinclude_HEADERS = NEUTRO.hxx + diff --git a/src/NEUTRO/NEUTRO.cxx b/src/NEUTRO/NEUTRO.cxx new file mode 100755 index 0000000..45494c0 --- /dev/null +++ b/src/NEUTRO/NEUTRO.cxx @@ -0,0 +1,157 @@ +using namespace std; +#include "NEUTRO.hxx" +#include +#include + +#include +#include + +extern "C" void transit_(void *compo,double *dt); +extern "C" void perma_(void *compo); + +//! Constructor for component "NEUTRO" instance +/*! + * + */ +NEUTRO_i::NEUTRO_i(CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + : Superv_Component_i(orb, poa, contId, instanceName, interfaceName) +{ + cerr << "create component" << endl; + _thisObj = this ; + _id = _poa->activate_object(_thisObj); +} + +//! Destructor for component "NEUTRO" instance +NEUTRO_i::~NEUTRO_i() +{ +} + +//! Register datastream ports for a component service given its name +/*! + * \param service_name : service name + * \return true if port registering succeeded, false if not + */ +CORBA::Boolean +NEUTRO_i::init_service(const char * service_name) { + CORBA::Boolean rtn = false; + string s_name(service_name); + if (s_name == "prun") + { + try + { + std::cerr << "NEUTRO: prun: " << std::endl; + //initialization CALCIUM ports IN + calcium_real_port_provides * p1; + p1 = add_port("CALCIUM_real","provides","tempi"); + p1->setDependencyType(CalciumTypes::ITERATION_DEPENDENCY); + calcium_integer_port_provides * p2; + p2 = add_port("CALCIUM_integer","provides","iconv"); + p2->setDependencyType(CalciumTypes::ITERATION_DEPENDENCY); + //initialization CALCIUM ports OUT + add_port("CALCIUM_real","uses","puissi"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "NEUTRO: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "NEUTRO: unknown exception" << std::endl; + } + + rtn = true; + } + if (s_name == "trun") + { + try + { + std::cerr << "NEUTRO: trun: " << std::endl; + //initialization CALCIUM ports IN + calcium_real_port_provides * p1; + p1 = add_port("CALCIUM_real","provides","temperature"); + p1->setDependencyType(CalciumTypes::TIME_DEPENDENCY); + //initialization CALCIUM ports OUT + add_port("CALCIUM_real","uses","puissance"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "NEUTRO: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "NEUTRO: unknown exception" << std::endl; + } + rtn = true; + } + return rtn; +} + +void NEUTRO_i::prun() +{ + std::cerr << "NEUTRO_i::prun" << std::endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + try + { + perma_(&component); + //to do or not to do ??? + cp_fin(component,CP_ARRET); + } + catch ( const CalciumException & ex) + { + std::cerr << ex.what() << std::endl; + cp_fin(component,CP_ARRET); + } + catch (...) + { + std::cerr << "unexpected exception" << std::endl; + cp_fin(component,CP_ARRET); + } + std::cerr << "end of NEUTRO_i::prun" << std::endl; +} + +void NEUTRO_i::trun(CORBA::Double dt) +{ + std::cerr << "NEUTRO_i::trun" << std::endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + try + { + transit_(&component,&dt); + //to do or not to do ??? + cp_fin(component,CP_ARRET); + } + catch ( const CalciumException & ex) + { + std::cerr << ex.what() << std::endl; + cp_fin(component,CP_ARRET); + } + catch (...) + { + std::cerr << "unexpected exception" << std::endl; + cp_fin(component,CP_ARRET); + } + std::cerr << "end of NEUTRO_i::trun" << std::endl; +} + +extern "C" +{ + PortableServer::ObjectId * NEUTROEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + MESSAGE("PortableServer::ObjectId * NEUTROEngine_factory()"); + NEUTRO_i * myEngine = new NEUTRO_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } +} diff --git a/src/NEUTRO/NEUTRO.hxx b/src/NEUTRO/NEUTRO.hxx new file mode 100644 index 0000000..94fea15 --- /dev/null +++ b/src/NEUTRO/NEUTRO.hxx @@ -0,0 +1,29 @@ +#ifndef _NEUTRO_HXX_ +#define _NEUTRO_HXX_ + +#include "Superv_Component_i.hxx" +#include "DSCCODE.hh" + +class NEUTRO_i: + public virtual POA_DSCCODE::NEUTRO, + public virtual Superv_Component_i +{ + public: + NEUTRO_i(CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, const char *interfaceName); + virtual ~NEUTRO_i(); + CORBA::Boolean init_service(const char * service_name); + void prun(); + void trun(CORBA::Double dt); +}; + +extern "C" +{ + PortableServer::ObjectId * NEUTROEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName); +} +#endif diff --git a/src/NEUTRO/neutro.f b/src/NEUTRO/neutro.f new file mode 100644 index 0000000..682f524 --- /dev/null +++ b/src/NEUTRO/neutro.f @@ -0,0 +1,162 @@ + subroutine transit(compo,dt0) + include "calcium.hf" + dimension q(24),temp(24) + double precision dt0 + integer compo + + dt=dt0 + + ti=0. + + do j=1,6 + npt=j*4 + q(npt)=0. + enddo + + pos=1. + do i=1,3 + do j=1,6 + npt=i+(j-1)*4 + haut=float(j)/6. + if(haut.gt.pos)then + q(npt)=0. + else + q(npt)=100. + endif + enddo + enddo +c +c initialisation puissance a t=0 +c + CALL cpeRE(compo,CP_TEMPS, ti, 1, 'puissance', 24, + & q , info) + IF( info.NE. CPOK )GO TO 9000 + + do while( .TRUE. ) +c do while( ti.lT.100. ) + + tf=ti+dt +c +c lecture de la temperature combustible entre ti et tf +c + CALL cplRE(compo,CP_TEMPS, ti,tf, npas, 'temperature', 24, + & nval, temp , info) + IF( info.NE. CPOK )GO TO 9000 + + do j=1,6 + npt=j*4 + q(npt)=0. + enddo +c +c calcul de la puissance degagee en fonction de la position +c des barres et de la temperature +c + do i=1,3 + do j=1,6 + npt=i+(j-1)*4 + haut=float(j)/6. + if(haut.gt.pos)then + q(npt)=0. + else + q(npt)=100.*(1.-0.0001*(temp(npt)-1000.)) + endif + enddo + enddo + write(6,*)"NEUTRO:","temps=",tf +c +c ecriture de la puissance a tf +c + CALL cpeRE(compo,CP_TEMPS, tf, npas+1, 'puissance', 24, + & q , info) + IF( info.NE. CPOK )GO TO 9000 + ti=tf + + enddo +9000 continue + CALL cpfin(compo,CP_ARRET, info) + end + + subroutine perma(compo) + include "calcium.hf" + dimension q(24),temp(24) + integer compo + + do j=1,6 + npt=j*4 + q(npt)=0. + enddo + + pos=1. + do i=1,3 + do j=1,6 + npt=i+(j-1)*4 + haut=float(j)/6. + if(haut.gt.pos)then + q(npt)=0. + else + q(npt)=100. + endif + enddo + enddo + + iter=0 + iconv=0 +c +c initialisation puissance a iter=0 +c + CALL cpeRE(compo,CP_ITERATION, ti, iter, 'puissi', 24, + & q , info) + IF( info.NE. CPOK )GO TO 9000 + + do while( iconv .eq. 0) +c +c lecture de la temperature combustible a iter +c + CALL cplRE(compo,CP_ITERATION, ti,tf, iter, 'tempi', 24, + & nval, temp , info) + IF( info.NE. CPOK )GO TO 9000 + + do j=1,6 + npt=j*4 + q(npt)=0. + enddo + +c +c calcul de la puissance degagee en fonction de la position +c des barres et de la temperature +c + do i=1,3 + do j=1,6 + npt=i+(j-1)*4 + haut=float(j)/6. + if(haut.gt.pos)then + q(npt)=0. + else + q(npt)=100.*(1.-0.0001*(temp(npt)-1000.)) + endif + enddo + enddo +c +c ecriture de la puissance a iter+1 +c + iter=iter+1 + CALL cpeRE(compo,CP_ITERATION, ti, iter, 'puissi', 24, + & q , info) + IF( info.NE. CPOK )GO TO 9000 +c +c lecture du flag de convergence iconv +c + CALL cplEN(compo,CP_ITERATION,ti, tf, iter , 'iconv', 1, + & nval, iconv , info) + write(6,*)"info:",info + write(6,*)"NEUTRO:",iter,iconv + CALL FLUSH(6) + IF( info.NE. CPOK )GO TO 9000 + + if(iconv.eq.1)go to 9000 + + enddo + +9000 continue + CALL cpfin(compo,CP_ARRET, info) + end diff --git a/src/SOLIDE/Makefile.am b/src/SOLIDE/Makefile.am new file mode 100644 index 0000000..25a8b7f --- /dev/null +++ b/src/SOLIDE/Makefile.am @@ -0,0 +1,12 @@ +include $(top_srcdir)/adm_local/make_common_starter.am + +AM_CFLAGS=$(KERNEL_INCLUDES) -fexceptions + +salomelib_LTLIBRARIES = libSOLIDEEngine.la +libSOLIDEEngine_la_SOURCES = SOLIDE.cxx solid.f +nodist_libSOLIDEEngine_la_SOURCES = +libSOLIDEEngine_la_CXXFLAGS = -I$(top_builddir)/idl $(KERNEL_INCLUDES) +libSOLIDEEngine_la_FFLAGS = $(KERNEL_INCLUDES) -fexceptions +libSOLIDEEngine_la_LIBADD = -L$(top_builddir)/idl -lDSCCODE -lg2c +salomeinclude_HEADERS = SOLIDE.hxx + diff --git a/src/SOLIDE/SOLIDE.cxx b/src/SOLIDE/SOLIDE.cxx new file mode 100755 index 0000000..5cc2051 --- /dev/null +++ b/src/SOLIDE/SOLIDE.cxx @@ -0,0 +1,164 @@ +using namespace std; +#include "SOLIDE.hxx" +#include +#include + +#include +#include + +extern "C" void transit_(void *compo,double *dt); +extern "C" void perma_(void *compo); + +//! Constructor for component "SOLIDE" instance +/*! + * + */ +SOLIDE_i::SOLIDE_i(CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + : Superv_Component_i(orb, poa, contId, instanceName, interfaceName) +{ + cerr << "create component" << endl; + _thisObj = this ; + _id = _poa->activate_object(_thisObj); +} + +//! Destructor for component "SOLIDE" instance +SOLIDE_i::~SOLIDE_i() +{ +} + +//! Register datastream ports for a component service given its name +/*! + * \param service_name : service name + * \return true if port registering succeeded, false if not + */ +CORBA::Boolean +SOLIDE_i::init_service(const char * service_name) { + CORBA::Boolean rtn = false; + string s_name(service_name); + if (s_name == "prun") + { + try + { + std::cerr << "SOLIDE: prun: " << std::endl; + //initialization CALCIUM ports IN + calcium_real_port_provides * p1; + p1 = add_port("CALCIUM_real","provides","puissi"); + p1->setDependencyType(CalciumTypes::ITERATION_DEPENDENCY); + p1 = add_port("CALCIUM_real","provides","tfi"); + p1->setDependencyType(CalciumTypes::ITERATION_DEPENDENCY); + //initialization CALCIUM ports OUT + add_port("CALCIUM_real","uses","tempi"); + add_port("CALCIUM_real","uses","tpi"); + add_port("CALCIUM_integer","uses","iconv"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "SOLIDE: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "SOLIDE: unknown exception" << std::endl; + } + + rtn = true; + } + if (s_name == "trun") + { + try + { + std::cerr << "SOLIDE: trun: " << std::endl; + //initialization CALCIUM ports IN + calcium_real_port_provides * p1; + p1 = add_port("CALCIUM_real","provides","rext"); + p1->setDependencyType(CalciumTypes::TIME_DEPENDENCY); + p1 = add_port("CALCIUM_real","provides","puissa"); + p1->setDependencyType(CalciumTypes::TIME_DEPENDENCY); + p1 = add_port("CALCIUM_real","provides","text"); + p1->setDependencyType(CalciumTypes::TIME_DEPENDENCY); + //initialization CALCIUM ports OUT + add_port("CALCIUM_real","uses","tparoi"); + add_port("CALCIUM_real","uses","rparoi"); + add_port("CALCIUM_real","uses","temperature"); + } + catch(const PortAlreadyDefined& ex) + { + std::cerr << "SOLIDE: " << ex.what() << std::endl; + //Ports already created : we use them + } + catch ( ... ) + { + std::cerr << "SOLIDE: unknown exception" << std::endl; + } + rtn = true; + } + return rtn; +} + +void SOLIDE_i::prun() +{ + std::cerr << "SOLIDE_i::prun" << std::endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + try + { + perma_(&component); + //to do or not to do ??? + cp_fin(component,CP_ARRET); + } + catch ( const CalciumException & ex) + { + std::cerr << ex.what() << std::endl; + cp_fin(component,CP_ARRET); + } + catch (...) + { + std::cerr << "unexpected exception" << std::endl; + cp_fin(component,CP_ARRET); + } + std::cerr << "end of SOLIDE_i::prun" << std::endl; +} + +void SOLIDE_i::trun(CORBA::Double dt) +{ + std::cerr << "SOLIDE_i::trun" << std::endl; + Superv_Component_i * component = dynamic_cast(this); + char nom_instance[INSTANCE_LEN]; + int info = cp_cd(component,nom_instance); + try + { + transit_(&component,&dt); + //to do or not to do ??? + cp_fin(component,CP_ARRET); + } + catch ( const CalciumException & ex) + { + std::cerr << ex.what() << std::endl; + cp_fin(component,CP_ARRET); + } + catch (...) + { + std::cerr << "unexpected exception" << std::endl; + cp_fin(component,CP_ARRET); + } + std::cerr << "end of SOLIDE_i::trun" << std::endl; +} + +extern "C" +{ + PortableServer::ObjectId * SOLIDEEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName) + { + MESSAGE("PortableServer::ObjectId * SOLIDEEngine_factory()"); + SOLIDE_i * myEngine = new SOLIDE_i(orb, poa, contId, instanceName, interfaceName); + return myEngine->getId() ; + } +} diff --git a/src/SOLIDE/SOLIDE.hxx b/src/SOLIDE/SOLIDE.hxx new file mode 100644 index 0000000..509a2e2 --- /dev/null +++ b/src/SOLIDE/SOLIDE.hxx @@ -0,0 +1,29 @@ +#ifndef _SOLIDE_HXX_ +#define _SOLIDE_HXX_ + +#include "Superv_Component_i.hxx" +#include "DSCCODE.hh" + +class SOLIDE_i: + public virtual POA_DSCCODE::SOLIDE, + public virtual Superv_Component_i +{ + public: + SOLIDE_i(CORBA::ORB_ptr orb, PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, const char *interfaceName); + virtual ~SOLIDE_i(); + CORBA::Boolean init_service(const char * service_name); + void prun(); + void trun(CORBA::Double dt); +}; + +extern "C" +{ + PortableServer::ObjectId * SOLIDEEngine_factory( CORBA::ORB_ptr orb, + PortableServer::POA_ptr poa, + PortableServer::ObjectId * contId, + const char *instanceName, + const char *interfaceName); +} +#endif diff --git a/src/SOLIDE/solid.f b/src/SOLIDE/solid.f new file mode 100644 index 0000000..c1a00fd --- /dev/null +++ b/src/SOLIDE/solid.f @@ -0,0 +1,3801 @@ + subroutine transit(compo,dt0) + include "calcium.hf" + double precision a(5,24),b(24),tn(24) + double precision q(6) + dimension maille(2),t(24),puiss(24) + dimension text(6),rflu(6) + dimension tpar(6),rpar(6) + double precision dt0 + integer compo + + dt=dt0 + + ldab=5 + ldb=24 + nrhs=1 + npt=0 + r=1. + ro=1. + rext=0.5 + te=10. + do i=1,24 + tn(i)=100. + t(i)=100. + puiss(i)=100. + enddo + do j=1,6 + q(j)=50. + enddo +c +c construction de la matrice (laplacien) +c + do i = 1, 5 + do j = 1, 24 + a(i,j) = 0. + enddo + enddo + + do i=2,3 + do j=2,5 + npt=i+(j-1)*4 + a(1,npt)=4./r+ro/dt + a(2,npt)=-1./r + a(5,npt)=-1./r + enddo + enddo + do i=2,3 + npt=i + a(1,npt)=3./r+ro/dt + a(2,npt)=-1./r + a(5,npt)=-1./r + npt=i+20 + a(1,npt)=3./r+ro/dt + a(2,npt)=-1./r + a(5,npt)=0. + enddo + do j=2,5 + npt=1+(j-1)*4 + a(1,npt)=3./r+ro/dt + a(2,npt)=-1./r + a(5,npt)=-1./r + npt=4+(j-1)*4 + a(1,npt)=3./r+1./(r/2.+rext)+ro/dt + a(2,npt)=0. + a(5,npt)=-1./r + enddo + i=1 + a(1,i)=2./r+ro/dt + a(2,i)=-1./r + a(5,i)=-1./r + i=21 + a(1,i)=2./r+ro/dt + a(2,i)=-1./r + a(5,i)=-1./r + i=24 + a(1,i)=2./r+1./(r/2.+rext)+ro/dt + a(2,i)=-1./r + a(5,i)=-1./r + i=4 + a(1,i)=2./r+1./(r/2.+rext)+ro/dt + a(2,i)=0. + a(5,i)=-1./r +c +c factorisation de la matrice +c + n=24 + kd=4 + + call dPBTRF( 'L' , N, KD, A, LDAB, INFO ) + + maille(1)=4 + maille(2)=6 + + ti=0. + + do i=1,6 + tpar(i)=t(4*i) + rpar(i)=r/2. + enddo +c +c initialisation de la temperature a t=0 +c + CALL cpeRE(compo,CP_TEMPS, ti, npas, 'temperature', 24, + & t , info) + IF( info.NE. CPOK )GO TO 9000 +c +c initialisation de la temperature de bord a t=0. +c + CALL cpeRE(compo,CP_TEMPS, ti, npas, 'tparoi', 6, + & tpar , info) + IF( info.NE. CPOK )GO TO 9000 +c +c initialisation de la resistance thermique de bord a t=0. +c + CALL cpeRE(compo,CP_TEMPS, ti, npas, 'rparoi', 6, + & rpar , info) + IF( info.NE. CPOK )GO TO 9000 +c +c boucle temporelle infinie +c + do while( .TRUE. ) +c do while( ti.lT.100. ) + tf=ti+dt +c +c lecture de la puissance combustible entre ti et tf +c + CALL cplRE(compo,CP_TEMPS,ti, tf, npas, 'puissa', 24, + & nval, puiss , info) + IF( info.NE. CPOK )GO TO 9000 +c +c lecture de la temperature exterieure entre ti et tf +c + CALL cplRE(compo,CP_TEMPS,ti, tf, npas, 'text', 6, + & nval, text , info) + IF( info.NE. CPOK )GO TO 9000 +c +c lecture de la resistance exterieure entre ti et tf +c + CALL cplRE(compo,CP_TEMPS,ti, tf, npas, 'rext', 6, + & nval, rflu , info) + IF( info.NE. CPOK )GO TO 9000 +c +c calcul du second membre +c + do i=1,24 + b(i)=ro/dt*tn(i) + enddo + + do j=1,6 + npt=j*4 + b(npt)=b(npt)+text(j)/(r/2+rflu(j)) + enddo + + do npt=1,24 + b(npt)=b(npt)+puiss(npt) + enddo +c +c resolution du systeme lineaire +c + call dPBTRs( 'L' , N,kd, nrhs,A,LDAB,b,ldb,INFO ) + + do i=1,24 + tn(i)=b(i) + t(i)=b(i) + enddo + + do i=1,6 + tpar(i)=t(4*i) + rpar(i)=r/2. + enddo + write(6,*)"SOLIDE: temps=",tf + call flush(6) +c +c ecriture de la temperature a t=tf +c + CALL cpeRE(compo,CP_TEMPS, tf, npas, 'temperature', 24, + & t , info) + IF( info.NE. CPOK )GO TO 9000 +c +c ecriture de la temperature de paroi a t=tf +c + CALL cpeRE(compo,CP_TEMPS, tf, npas, 'tparoi', 6, + & tpar , info) + IF( info.NE. CPOK )GO TO 9000 +c +c ecriture de la resistance de bord a t=tf +c + CALL cpeRE(compo,CP_TEMPS, tf, npas, 'rparoi', 6, + & rpar , info) + IF( info.NE. CPOK )GO TO 9000 + + ti=tf + + enddo +9000 continue + CALL cpfin(compo,CP_ARRET, info) + end +c + SUBROUTINE DPBTRF( UPLO, N, KD, AB, LDAB, INFO ) +* +* -- LAPACK routine (version 2.0) -- +* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., +* Courant Institute, Argonne National Lab, and Rice University +* March 31, 1993 +* +* .. Scalar Arguments .. + CHARACTER UPLO + INTEGER INFO, KD, LDAB, N +* .. +* .. Array Arguments .. + DOUBLE PRECISION AB( LDAB, * ) +* .. +* +* Purpose +* ======= +* +* DPBTRF computes the Cholesky factorization of a real symmetric +* positive definite band matrix A. +* +* The factorization has the form +* A = U**T * U, if UPLO = 'U', or +* A = L * L**T, if UPLO = 'L', +* where U is an upper triangular matrix and L is lower triangular. +* +* Arguments +* ========= +* +* UPLO (input) CHARACTER*1 +* = 'U': Upper triangle of A is stored; +* = 'L': Lower triangle of A is stored. +* +* N (input) INTEGER +* The order of the matrix A. N >= 0. +* +* KD (input) INTEGER +* The number of superdiagonals of the matrix A if UPLO = 'U', +* or the number of subdiagonals if UPLO = 'L'. KD >= 0. +* +* AB (input/output) DOUBLE PRECISION array, dimension (LDAB,N) +* On entry, the upper or lower triangle of the symmetric band +* matrix A, stored in the first KD+1 rows of the array. The +* j-th column of A is stored in the j-th column of the array AB +* as follows: +* if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j; +* if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd). +* +* On exit, if INFO = 0, the triangular factor U or L from the +* Cholesky factorization A = U**T*U or A = L*L**T of the band +* matrix A, in the same storage format as A. +* +* LDAB (input) INTEGER +* The leading dimension of the array AB. LDAB >= KD+1. +* +* INFO (output) INTEGER +* = 0: successful exit +* < 0: if INFO = -i, the i-th argument had an illegal value +* > 0: if INFO = i, the leading minor of order i is not +* positive definite, and the factorization could not be +* completed. +* +* Further Details +* =============== +* +* The band storage scheme is illustrated by the following example, when +* N = 6, KD = 2, and UPLO = 'U': +* +* On entry: On exit: +* +* * * a13 a24 a35 a46 * * u13 u24 u35 u46 +* * a12 a23 a34 a45 a56 * u12 u23 u34 u45 u56 +* a11 a22 a33 a44 a55 a66 u11 u22 u33 u44 u55 u66 +* +* Similarly, if UPLO = 'L' the format of A is as follows: +* +* On entry: On exit: +* +* a11 a22 a33 a44 a55 a66 l11 l22 l33 l44 l55 l66 +* a21 a32 a43 a54 a65 * l21 l32 l43 l54 l65 * +* a31 a42 a53 a64 * * l31 l42 l53 l64 * * +* +* Array elements marked * are not used by the routine. +* +* Contributed by +* Peter Mayes and Giuseppe Radicati, IBM ECSEC, Rome, March 23, 1989 +* +* ===================================================================== +* +* .. Parameters .. + DOUBLE PRECISION ONE, ZERO + PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) + INTEGER NBMAX, LDWORK + PARAMETER ( NBMAX = 32, LDWORK = NBMAX+1 ) +* .. +* .. Local Scalars .. + INTEGER I, I2, I3, IB, II, J, JJ, NB +* .. +* .. Local Arrays .. + DOUBLE PRECISION WORK( LDWORK, NBMAX ) +* .. +* .. External Functions .. + LOGICAL LSAME + INTEGER ILAENV + EXTERNAL LSAME, ILAENV +* .. +* .. External Subroutines .. + EXTERNAL DGEMM, DPBTF2, DPOTF2, DSYRK, DTRSM, XERBLA +* .. +* .. Intrinsic Functions .. + INTRINSIC MIN +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* + INFO = 0 + IF( ( .NOT.LSAME( UPLO, 'U' ) ) .AND. + $ ( .NOT.LSAME( UPLO, 'L' ) ) ) THEN + INFO = -1 + ELSE IF( N.LT.0 ) THEN + INFO = -2 + ELSE IF( KD.LT.0 ) THEN + INFO = -3 + ELSE IF( LDAB.LT.KD+1 ) THEN + INFO = -5 + END IF + IF( INFO.NE.0 ) THEN + CALL XERBLA( 'DPBTRF', -INFO ) + RETURN + END IF +* +* Quick return if possible +* + IF( N.EQ.0 ) + $ RETURN +* +* Determine the block size for this environment +* + NB = ILAENV( 1, 'DPBTRF', UPLO, N, KD, -1, -1 ) +* +* The block size must not exceed the semi-bandwidth KD, and must not +* exceed the limit set by the size of the local array WORK. +* + NB = MIN( NB, NBMAX ) +* + IF( NB.LE.1 .OR. NB.GT.KD ) THEN +* +* Use unblocked code +* + CALL DPBTF2( UPLO, N, KD, AB, LDAB, INFO ) + ELSE +* +* Use blocked code +* + IF( LSAME( UPLO, 'U' ) ) THEN +* +* Compute the Cholesky factorization of a symmetric band +* matrix, given the upper triangle of the matrix in band +* storage. +* +* Zero the upper triangle of the work array. +* + DO 20 J = 1, NB + DO 10 I = 1, J - 1 + WORK( I, J ) = ZERO + 10 CONTINUE + 20 CONTINUE +* +* Process the band matrix one diagonal block at a time. +* + DO 70 I = 1, N, NB + IB = MIN( NB, N-I+1 ) +* +* Factorize the diagonal block +* + CALL DPOTF2( UPLO, IB, AB( KD+1, I ), LDAB-1, II ) + IF( II.NE.0 ) THEN + INFO = I + II - 1 + GO TO 150 + END IF + IF( I+IB.LE.N ) THEN +* +* Update the relevant part of the trailing submatrix. +* If A11 denotes the diagonal block which has just been +* factorized, then we need to update the remaining +* blocks in the diagram: +* +* A11 A12 A13 +* A22 A23 +* A33 +* +* The numbers of rows and columns in the partitioning +* are IB, I2, I3 respectively. The blocks A12, A22 and +* A23 are empty if IB = KD. The upper triangle of A13 +* lies outside the band. +* + I2 = MIN( KD-IB, N-I-IB+1 ) + I3 = MIN( IB, N-I-KD+1 ) +* + IF( I2.GT.0 ) THEN +* +* Update A12 +* + CALL DTRSM( 'Left', 'Upper', 'Transpose', + $ 'Non-unit', IB, I2, ONE, AB( KD+1, I ), + $ LDAB-1, AB( KD+1-IB, I+IB ), LDAB-1 ) +* +* Update A22 +* + CALL DSYRK( 'Upper', 'Transpose', I2, IB, -ONE, + $ AB( KD+1-IB, I+IB ), LDAB-1, ONE, + $ AB( KD+1, I+IB ), LDAB-1 ) + END IF +* + IF( I3.GT.0 ) THEN +* +* Copy the lower triangle of A13 into the work array. +* + DO 40 JJ = 1, I3 + DO 30 II = JJ, IB + WORK( II, JJ ) = AB( II-JJ+1, JJ+I+KD-1 ) + 30 CONTINUE + 40 CONTINUE +* +* Update A13 (in the work array). +* + CALL DTRSM( 'Left', 'Upper', 'Transpose', + $ 'Non-unit', IB, I3, ONE, AB( KD+1, I ), + $ LDAB-1, WORK, LDWORK ) +* +* Update A23 +* + IF( I2.GT.0 ) + $ CALL DGEMM( 'Transpose', 'No Transpose', I2, I3, + $ IB, -ONE, AB( KD+1-IB, I+IB ), + $ LDAB-1, WORK, LDWORK, ONE, + $ AB( 1+IB, I+KD ), LDAB-1 ) +* +* Update A33 +* + CALL DSYRK( 'Upper', 'Transpose', I3, IB, -ONE, + $ WORK, LDWORK, ONE, AB( KD+1, I+KD ), + $ LDAB-1 ) +* +* Copy the lower triangle of A13 back into place. +* + DO 60 JJ = 1, I3 + DO 50 II = JJ, IB + AB( II-JJ+1, JJ+I+KD-1 ) = WORK( II, JJ ) + 50 CONTINUE + 60 CONTINUE + END IF + END IF + 70 CONTINUE + ELSE +* +* Compute the Cholesky factorization of a symmetric band +* matrix, given the lower triangle of the matrix in band +* storage. +* +* Zero the lower triangle of the work array. +* + DO 90 J = 1, NB + DO 80 I = J + 1, NB + WORK( I, J ) = ZERO + 80 CONTINUE + 90 CONTINUE +* +* Process the band matrix one diagonal block at a time. +* + DO 140 I = 1, N, NB + IB = MIN( NB, N-I+1 ) +* +* Factorize the diagonal block +* + CALL DPOTF2( UPLO, IB, AB( 1, I ), LDAB-1, II ) + IF( II.NE.0 ) THEN + INFO = I + II - 1 + GO TO 150 + END IF + IF( I+IB.LE.N ) THEN +* +* Update the relevant part of the trailing submatrix. +* If A11 denotes the diagonal block which has just been +* factorized, then we need to update the remaining +* blocks in the diagram: +* +* A11 +* A21 A22 +* A31 A32 A33 +* +* The numbers of rows and columns in the partitioning +* are IB, I2, I3 respectively. The blocks A21, A22 and +* A32 are empty if IB = KD. The lower triangle of A31 +* lies outside the band. +* + I2 = MIN( KD-IB, N-I-IB+1 ) + I3 = MIN( IB, N-I-KD+1 ) +* + IF( I2.GT.0 ) THEN +* +* Update A21 +* + CALL DTRSM( 'Right', 'Lower', 'Transpose', + $ 'Non-unit', I2, IB, ONE, AB( 1, I ), + $ LDAB-1, AB( 1+IB, I ), LDAB-1 ) +* +* Update A22 +* + CALL DSYRK( 'Lower', 'No Transpose', I2, IB, -ONE, + $ AB( 1+IB, I ), LDAB-1, ONE, + $ AB( 1, I+IB ), LDAB-1 ) + END IF +* + IF( I3.GT.0 ) THEN +* +* Copy the upper triangle of A31 into the work array. +* + DO 110 JJ = 1, IB + DO 100 II = 1, MIN( JJ, I3 ) + WORK( II, JJ ) = AB( KD+1-JJ+II, JJ+I-1 ) + 100 CONTINUE + 110 CONTINUE +* +* Update A31 (in the work array). +* + CALL DTRSM( 'Right', 'Lower', 'Transpose', + $ 'Non-unit', I3, IB, ONE, AB( 1, I ), + $ LDAB-1, WORK, LDWORK ) +* +* Update A32 +* + IF( I2.GT.0 ) + $ CALL DGEMM( 'No transpose', 'Transpose', I3, I2, + $ IB, -ONE, WORK, LDWORK, + $ AB( 1+IB, I ), LDAB-1, ONE, + $ AB( 1+KD-IB, I+IB ), LDAB-1 ) +* +* Update A33 +* + CALL DSYRK( 'Lower', 'No Transpose', I3, IB, -ONE, + $ WORK, LDWORK, ONE, AB( 1, I+KD ), + $ LDAB-1 ) +* +* Copy the upper triangle of A31 back into place. +* + DO 130 JJ = 1, IB + DO 120 II = 1, MIN( JJ, I3 ) + AB( KD+1-JJ+II, JJ+I-1 ) = WORK( II, JJ ) + 120 CONTINUE + 130 CONTINUE + END IF + END IF + 140 CONTINUE + END IF + END IF + RETURN +* + 150 CONTINUE + RETURN +* +* End of DPBTRF +* + END + SUBROUTINE DPBTRS( UPLO, N, KD, NRHS, AB, LDAB, B, LDB, INFO ) +* +* -- LAPACK routine (version 2.0) -- +* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., +* Courant Institute, Argonne National Lab, and Rice University +* September 30, 1994 +* +* .. Scalar Arguments .. + CHARACTER UPLO + INTEGER INFO, KD, LDAB, LDB, N, NRHS +* .. +* .. Array Arguments .. + DOUBLE PRECISION AB( LDAB, * ), B( LDB, * ) +* .. +* +* Purpose +* ======= +* +* DPBTRS solves a system of linear equations A*X = B with a symmetric +* positive definite band matrix A using the Cholesky factorization +* A = U**T*U or A = L*L**T computed by DPBTRF. +* +* Arguments +* ========= +* +* UPLO (input) CHARACTER*1 +* = 'U': Upper triangular factor stored in AB; +* = 'L': Lower triangular factor stored in AB. +* +* N (input) INTEGER +* The order of the matrix A. N >= 0. +* +* KD (input) INTEGER +* The number of superdiagonals of the matrix A if UPLO = 'U', +* or the number of subdiagonals if UPLO = 'L'. KD >= 0. +* +* NRHS (input) INTEGER +* The number of right hand sides, i.e., the number of columns +* of the matrix B. NRHS >= 0. +* +* AB (input) DOUBLE PRECISION array, dimension (LDAB,N) +* The triangular factor U or L from the Cholesky factorization +* A = U**T*U or A = L*L**T of the band matrix A, stored in the +* first KD+1 rows of the array. The j-th column of U or L is +* stored in the j-th column of the array AB as follows: +* if UPLO ='U', AB(kd+1+i-j,j) = U(i,j) for max(1,j-kd)<=i<=j; +* if UPLO ='L', AB(1+i-j,j) = L(i,j) for j<=i<=min(n,j+kd). +* +* LDAB (input) INTEGER +* The leading dimension of the array AB. LDAB >= KD+1. +* +* B (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS) +* On entry, the right hand side matrix B. +* On exit, the solution matrix X. +* +* LDB (input) INTEGER +* The leading dimension of the array B. LDB >= max(1,N). +* +* INFO (output) INTEGER +* = 0: successful exit +* < 0: if INFO = -i, the i-th argument had an illegal value +* +* ===================================================================== +* +* .. Local Scalars .. + LOGICAL UPPER + INTEGER J +* .. +* .. External Functions .. + LOGICAL LSAME + EXTERNAL LSAME +* .. +* .. External Subroutines .. + EXTERNAL DTBSV, XERBLA +* .. +* .. Intrinsic Functions .. + INTRINSIC MAX +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* + INFO = 0 + UPPER = LSAME( UPLO, 'U' ) + IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN + INFO = -1 + ELSE IF( N.LT.0 ) THEN + INFO = -2 + ELSE IF( KD.LT.0 ) THEN + INFO = -3 + ELSE IF( NRHS.LT.0 ) THEN + INFO = -4 + ELSE IF( LDAB.LT.KD+1 ) THEN + INFO = -6 + ELSE IF( LDB.LT.MAX( 1, N ) ) THEN + INFO = -8 + END IF + IF( INFO.NE.0 ) THEN + CALL XERBLA( 'DPBTRS', -INFO ) + RETURN + END IF +* +* Quick return if possible +* + IF( N.EQ.0 .OR. NRHS.EQ.0 ) + $ RETURN +* + IF( UPPER ) THEN +* +* Solve A*X = B where A = U'*U. +* + DO 10 J = 1, NRHS +* +* Solve U'*X = B, overwriting B with X. +* + CALL DTBSV( 'Upper', 'Transpose', 'Non-unit', N, KD, AB, + $ LDAB, B( 1, J ), 1 ) +* +* Solve U*X = B, overwriting B with X. +* + CALL DTBSV( 'Upper', 'No transpose', 'Non-unit', N, KD, AB, + $ LDAB, B( 1, J ), 1 ) + 10 CONTINUE + ELSE +* +* Solve A*X = B where A = L*L'. +* + DO 20 J = 1, NRHS +* +* Solve L*X = B, overwriting B with X. +* + CALL DTBSV( 'Lower', 'No transpose', 'Non-unit', N, KD, AB, + $ LDAB, B( 1, J ), 1 ) +* +* Solve L'*X = B, overwriting B with X. +* + CALL DTBSV( 'Lower', 'Transpose', 'Non-unit', N, KD, AB, + $ LDAB, B( 1, J ), 1 ) + 20 CONTINUE + END IF +* + RETURN +* +* End of DPBTRS +* + END + INTEGER FUNCTION ILAENV( ISPEC, NAME, OPTS, N1, N2, N3, + $ N4 ) +* +* -- LAPACK auxiliary routine (version 2.0) -- +* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., +* Courant Institute, Argonne National Lab, and Rice University +* September 30, 1994 +* +* .. Scalar Arguments .. + CHARACTER*( * ) NAME, OPTS + INTEGER ISPEC, N1, N2, N3, N4 +* .. +* +* Purpose +* ======= +* +* ILAENV is called from the LAPACK routines to choose problem-dependent +* parameters for the local environment. See ISPEC for a description of +* the parameters. +* +* This version provides a set of parameters which should give good, +* but not optimal, performance on many of the currently available +* computers. Users are encouraged to modify this subroutine to set +* the tuning parameters for their particular machine using the option +* and problem size information in the arguments. +* +* This routine will not function correctly if it is converted to all +* lower case. Converting it to all upper case is allowed. +* +* Arguments +* ========= +* +* ISPEC (input) INTEGER +* Specifies the parameter to be returned as the value of +* ILAENV. +* = 1: the optimal blocksize; if this value is 1, an unblocked +* algorithm will give the best performance. +* = 2: the minimum block size for which the block routine +* should be used; if the usable block size is less than +* this value, an unblocked routine should be used. +* = 3: the crossover point (in a block routine, for N less +* than this value, an unblocked routine should be used) +* = 4: the number of shifts, used in the nonsymmetric +* eigenvalue routines +* = 5: the minimum column dimension for blocking to be used; +* rectangular blocks must have dimension at least k by m, +* where k is given by ILAENV(2,...) and m by ILAENV(5,...) +* = 6: the crossover point for the SVD (when reducing an m by n +* matrix to bidiagonal form, if max(m,n)/min(m,n) exceeds +* this value, a QR factorization is used first to reduce +* the matrix to a triangular form.) +* = 7: the number of processors +* = 8: the crossover point for the multishift QR and QZ methods +* for nonsymmetric eigenvalue problems. +* +* NAME (input) CHARACTER*(*) +* The name of the calling subroutine, in either upper case or +* lower case. +* +* OPTS (input) CHARACTER*(*) +* The character options to the subroutine NAME, concatenated +* into a single character string. For example, UPLO = 'U', +* TRANS = 'T', and DIAG = 'N' for a triangular routine would +* be specified as OPTS = 'UTN'. +* +* N1 (input) INTEGER +* N2 (input) INTEGER +* N3 (input) INTEGER +* N4 (input) INTEGER +* Problem dimensions for the subroutine NAME; these may not all +* be required. +* +* (ILAENV) (output) INTEGER +* >= 0: the value of the parameter specified by ISPEC +* < 0: if ILAENV = -k, the k-th argument had an illegal value. +* +* Further Details +* =============== +* +* The following conventions have been used when calling ILAENV from the +* LAPACK routines: +* 1) OPTS is a concatenation of all of the character options to +* subroutine NAME, in the same order that they appear in the +* argument list for NAME, even if they are not used in determining +* the value of the parameter specified by ISPEC. +* 2) The problem dimensions N1, N2, N3, N4 are specified in the order +* that they appear in the argument list for NAME. N1 is used +* first, N2 second, and so on, and unused problem dimensions are +* passed a value of -1. +* 3) The parameter value returned by ILAENV is checked for validity in +* the calling subroutine. For example, ILAENV is used to retrieve +* the optimal blocksize for STRTRI as follows: +* +* NB = ILAENV( 1, 'STRTRI', UPLO // DIAG, N, -1, -1, -1 ) +* IF( NB.LE.1 ) NB = MAX( 1, N ) +* +* ===================================================================== +* +* .. Local Scalars .. + LOGICAL CNAME, SNAME + CHARACTER*1 C1 + CHARACTER*2 C2, C4 + CHARACTER*3 C3 + CHARACTER*6 SUBNAM + INTEGER I, IC, IZ, NB, NBMIN, NX +* .. +* .. Intrinsic Functions .. + INTRINSIC CHAR, ICHAR, INT, MIN, REAL +* .. +* .. Executable Statements .. +* + GO TO ( 100, 100, 100, 400, 500, 600, 700, 800 ) ISPEC +* +* Invalid value for ISPEC +* + ILAENV = -1 + RETURN +* + 100 CONTINUE +* +* Convert NAME to upper case if the first character is lower case. +* + ILAENV = 1 + SUBNAM = NAME + IC = ICHAR( SUBNAM( 1:1 ) ) + IZ = ICHAR( 'Z' ) + IF( IZ.EQ.90 .OR. IZ.EQ.122 ) THEN +* +* ASCII character set +* + IF( IC.GE.97 .AND. IC.LE.122 ) THEN + SUBNAM( 1:1 ) = CHAR( IC-32 ) + DO 10 I = 2, 6 + IC = ICHAR( SUBNAM( I:I ) ) + IF( IC.GE.97 .AND. IC.LE.122 ) + $ SUBNAM( I:I ) = CHAR( IC-32 ) + 10 CONTINUE + END IF +* + ELSE IF( IZ.EQ.233 .OR. IZ.EQ.169 ) THEN +* +* EBCDIC character set +* + IF( ( IC.GE.129 .AND. IC.LE.137 ) .OR. + $ ( IC.GE.145 .AND. IC.LE.153 ) .OR. + $ ( IC.GE.162 .AND. IC.LE.169 ) ) THEN + SUBNAM( 1:1 ) = CHAR( IC+64 ) + DO 20 I = 2, 6 + IC = ICHAR( SUBNAM( I:I ) ) + IF( ( IC.GE.129 .AND. IC.LE.137 ) .OR. + $ ( IC.GE.145 .AND. IC.LE.153 ) .OR. + $ ( IC.GE.162 .AND. IC.LE.169 ) ) + $ SUBNAM( I:I ) = CHAR( IC+64 ) + 20 CONTINUE + END IF +* + ELSE IF( IZ.EQ.218 .OR. IZ.EQ.250 ) THEN +* +* Prime machines: ASCII+128 +* + IF( IC.GE.225 .AND. IC.LE.250 ) THEN + SUBNAM( 1:1 ) = CHAR( IC-32 ) + DO 30 I = 2, 6 + IC = ICHAR( SUBNAM( I:I ) ) + IF( IC.GE.225 .AND. IC.LE.250 ) + $ SUBNAM( I:I ) = CHAR( IC-32 ) + 30 CONTINUE + END IF + END IF +* + C1 = SUBNAM( 1:1 ) + SNAME = C1.EQ.'S' .OR. C1.EQ.'D' + CNAME = C1.EQ.'C' .OR. C1.EQ.'Z' + IF( .NOT.( CNAME .OR. SNAME ) ) + $ RETURN + C2 = SUBNAM( 2:3 ) + C3 = SUBNAM( 4:6 ) + C4 = C3( 2:3 ) +* + GO TO ( 110, 200, 300 ) ISPEC +* + 110 CONTINUE +* +* ISPEC = 1: block size +* +* In these examples, separate code is provided for setting NB for +* real and complex. We assume that NB will take the same value in +* single or double precision. +* + NB = 1 +* + IF( C2.EQ.'GE' ) THEN + IF( C3.EQ.'TRF' ) THEN + IF( SNAME ) THEN + NB = 64 + ELSE + NB = 64 + END IF + ELSE IF( C3.EQ.'QRF' .OR. C3.EQ.'RQF' .OR. C3.EQ.'LQF' .OR. + $ C3.EQ.'QLF' ) THEN + IF( SNAME ) THEN + NB = 32 + ELSE + NB = 32 + END IF + ELSE IF( C3.EQ.'HRD' ) THEN + IF( SNAME ) THEN + NB = 32 + ELSE + NB = 32 + END IF + ELSE IF( C3.EQ.'BRD' ) THEN + IF( SNAME ) THEN + NB = 32 + ELSE + NB = 32 + END IF + ELSE IF( C3.EQ.'TRI' ) THEN + IF( SNAME ) THEN + NB = 64 + ELSE + NB = 64 + END IF + END IF + ELSE IF( C2.EQ.'PO' ) THEN + IF( C3.EQ.'TRF' ) THEN + IF( SNAME ) THEN + NB = 64 + ELSE + NB = 64 + END IF + END IF + ELSE IF( C2.EQ.'SY' ) THEN + IF( C3.EQ.'TRF' ) THEN + IF( SNAME ) THEN + NB = 64 + ELSE + NB = 64 + END IF + ELSE IF( SNAME .AND. C3.EQ.'TRD' ) THEN + NB = 1 + ELSE IF( SNAME .AND. C3.EQ.'GST' ) THEN + NB = 64 + END IF + ELSE IF( CNAME .AND. C2.EQ.'HE' ) THEN + IF( C3.EQ.'TRF' ) THEN + NB = 64 + ELSE IF( C3.EQ.'TRD' ) THEN + NB = 1 + ELSE IF( C3.EQ.'GST' ) THEN + NB = 64 + END IF + ELSE IF( SNAME .AND. C2.EQ.'OR' ) THEN + IF( C3( 1:1 ).EQ.'G' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NB = 32 + END IF + ELSE IF( C3( 1:1 ).EQ.'M' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NB = 32 + END IF + END IF + ELSE IF( CNAME .AND. C2.EQ.'UN' ) THEN + IF( C3( 1:1 ).EQ.'G' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NB = 32 + END IF + ELSE IF( C3( 1:1 ).EQ.'M' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NB = 32 + END IF + END IF + ELSE IF( C2.EQ.'GB' ) THEN + IF( C3.EQ.'TRF' ) THEN + IF( SNAME ) THEN + IF( N4.LE.64 ) THEN + NB = 1 + ELSE + NB = 32 + END IF + ELSE + IF( N4.LE.64 ) THEN + NB = 1 + ELSE + NB = 32 + END IF + END IF + END IF + ELSE IF( C2.EQ.'PB' ) THEN + IF( C3.EQ.'TRF' ) THEN + IF( SNAME ) THEN + IF( N2.LE.64 ) THEN + NB = 1 + ELSE + NB = 32 + END IF + ELSE + IF( N2.LE.64 ) THEN + NB = 1 + ELSE + NB = 32 + END IF + END IF + END IF + ELSE IF( C2.EQ.'TR' ) THEN + IF( C3.EQ.'TRI' ) THEN + IF( SNAME ) THEN + NB = 64 + ELSE + NB = 64 + END IF + END IF + ELSE IF( C2.EQ.'LA' ) THEN + IF( C3.EQ.'UUM' ) THEN + IF( SNAME ) THEN + NB = 64 + ELSE + NB = 64 + END IF + END IF + ELSE IF( SNAME .AND. C2.EQ.'ST' ) THEN + IF( C3.EQ.'EBZ' ) THEN + NB = 1 + END IF + END IF + ILAENV = NB + RETURN +* + 200 CONTINUE +* +* ISPEC = 2: minimum block size +* + NBMIN = 2 + IF( C2.EQ.'GE' ) THEN + IF( C3.EQ.'QRF' .OR. C3.EQ.'RQF' .OR. C3.EQ.'LQF' .OR. + $ C3.EQ.'QLF' ) THEN + IF( SNAME ) THEN + NBMIN = 2 + ELSE + NBMIN = 2 + END IF + ELSE IF( C3.EQ.'HRD' ) THEN + IF( SNAME ) THEN + NBMIN = 2 + ELSE + NBMIN = 2 + END IF + ELSE IF( C3.EQ.'BRD' ) THEN + IF( SNAME ) THEN + NBMIN = 2 + ELSE + NBMIN = 2 + END IF + ELSE IF( C3.EQ.'TRI' ) THEN + IF( SNAME ) THEN + NBMIN = 2 + ELSE + NBMIN = 2 + END IF + END IF + ELSE IF( C2.EQ.'SY' ) THEN + IF( C3.EQ.'TRF' ) THEN + IF( SNAME ) THEN + NBMIN = 8 + ELSE + NBMIN = 8 + END IF + ELSE IF( SNAME .AND. C3.EQ.'TRD' ) THEN + NBMIN = 2 + END IF + ELSE IF( CNAME .AND. C2.EQ.'HE' ) THEN + IF( C3.EQ.'TRD' ) THEN + NBMIN = 2 + END IF + ELSE IF( SNAME .AND. C2.EQ.'OR' ) THEN + IF( C3( 1:1 ).EQ.'G' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NBMIN = 2 + END IF + ELSE IF( C3( 1:1 ).EQ.'M' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NBMIN = 2 + END IF + END IF + ELSE IF( CNAME .AND. C2.EQ.'UN' ) THEN + IF( C3( 1:1 ).EQ.'G' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NBMIN = 2 + END IF + ELSE IF( C3( 1:1 ).EQ.'M' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NBMIN = 2 + END IF + END IF + END IF + ILAENV = NBMIN + RETURN +* + 300 CONTINUE +* +* ISPEC = 3: crossover point +* + NX = 0 + IF( C2.EQ.'GE' ) THEN + IF( C3.EQ.'QRF' .OR. C3.EQ.'RQF' .OR. C3.EQ.'LQF' .OR. + $ C3.EQ.'QLF' ) THEN + IF( SNAME ) THEN + NX = 128 + ELSE + NX = 128 + END IF + ELSE IF( C3.EQ.'HRD' ) THEN + IF( SNAME ) THEN + NX = 128 + ELSE + NX = 128 + END IF + ELSE IF( C3.EQ.'BRD' ) THEN + IF( SNAME ) THEN + NX = 128 + ELSE + NX = 128 + END IF + END IF + ELSE IF( C2.EQ.'SY' ) THEN + IF( SNAME .AND. C3.EQ.'TRD' ) THEN + NX = 1 + END IF + ELSE IF( CNAME .AND. C2.EQ.'HE' ) THEN + IF( C3.EQ.'TRD' ) THEN + NX = 1 + END IF + ELSE IF( SNAME .AND. C2.EQ.'OR' ) THEN + IF( C3( 1:1 ).EQ.'G' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NX = 128 + END IF + END IF + ELSE IF( CNAME .AND. C2.EQ.'UN' ) THEN + IF( C3( 1:1 ).EQ.'G' ) THEN + IF( C4.EQ.'QR' .OR. C4.EQ.'RQ' .OR. C4.EQ.'LQ' .OR. + $ C4.EQ.'QL' .OR. C4.EQ.'HR' .OR. C4.EQ.'TR' .OR. + $ C4.EQ.'BR' ) THEN + NX = 128 + END IF + END IF + END IF + ILAENV = NX + RETURN +* + 400 CONTINUE +* +* ISPEC = 4: number of shifts (used by xHSEQR) +* + ILAENV = 6 + RETURN +* + 500 CONTINUE +* +* ISPEC = 5: minimum column dimension (not used) +* + ILAENV = 2 + RETURN +* + 600 CONTINUE +* +* ISPEC = 6: crossover point for SVD (used by xGELSS and xGESVD) +* + ILAENV = INT( REAL( MIN( N1, N2 ) )*1.6E0 ) + RETURN +* + 700 CONTINUE +* +* ISPEC = 7: number of processors (not used) +* + ILAENV = 1 + RETURN +* + 800 CONTINUE +* +* ISPEC = 8: crossover point for multishift (used by xHSEQR) +* + ILAENV = 50 + RETURN +* +* End of ILAENV +* + END + SUBROUTINE DPBTF2( UPLO, N, KD, AB, LDAB, INFO ) +* +* -- LAPACK routine (version 2.0) -- +* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., +* Courant Institute, Argonne National Lab, and Rice University +* February 29, 1992 +* +* .. Scalar Arguments .. + CHARACTER UPLO + INTEGER INFO, KD, LDAB, N +* .. +* .. Array Arguments .. + DOUBLE PRECISION AB( LDAB, * ) +* .. +* +* Purpose +* ======= +* +* DPBTF2 computes the Cholesky factorization of a real symmetric +* positive definite band matrix A. +* +* The factorization has the form +* A = U' * U , if UPLO = 'U', or +* A = L * L', if UPLO = 'L', +* where U is an upper triangular matrix, U' is the transpose of U, and +* L is lower triangular. +* +* This is the unblocked version of the algorithm, calling Level 2 BLAS. +* +* Arguments +* ========= +* +* UPLO (input) CHARACTER*1 +* Specifies whether the upper or lower triangular part of the +* symmetric matrix A is stored: +* = 'U': Upper triangular +* = 'L': Lower triangular +* +* N (input) INTEGER +* The order of the matrix A. N >= 0. +* +* KD (input) INTEGER +* The number of super-diagonals of the matrix A if UPLO = 'U', +* or the number of sub-diagonals if UPLO = 'L'. KD >= 0. +* +* AB (input/output) DOUBLE PRECISION array, dimension (LDAB,N) +* On entry, the upper or lower triangle of the symmetric band +* matrix A, stored in the first KD+1 rows of the array. The +* j-th column of A is stored in the j-th column of the array AB +* as follows: +* if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j; +* if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd). +* +* On exit, if INFO = 0, the triangular factor U or L from the +* Cholesky factorization A = U'*U or A = L*L' of the band +* matrix A, in the same storage format as A. +* +* LDAB (input) INTEGER +* The leading dimension of the array AB. LDAB >= KD+1. +* +* INFO (output) INTEGER +* = 0: successful exit +* < 0: if INFO = -k, the k-th argument had an illegal value +* > 0: if INFO = k, the leading minor of order k is not +* positive definite, and the factorization could not be +* completed. +* +* Further Details +* =============== +* +* The band storage scheme is illustrated by the following example, when +* N = 6, KD = 2, and UPLO = 'U': +* +* On entry: On exit: +* +* * * a13 a24 a35 a46 * * u13 u24 u35 u46 +* * a12 a23 a34 a45 a56 * u12 u23 u34 u45 u56 +* a11 a22 a33 a44 a55 a66 u11 u22 u33 u44 u55 u66 +* +* Similarly, if UPLO = 'L' the format of A is as follows: +* +* On entry: On exit: +* +* a11 a22 a33 a44 a55 a66 l11 l22 l33 l44 l55 l66 +* a21 a32 a43 a54 a65 * l21 l32 l43 l54 l65 * +* a31 a42 a53 a64 * * l31 l42 l53 l64 * * +* +* Array elements marked * are not used by the routine. +* +* ===================================================================== +* +* .. Parameters .. + DOUBLE PRECISION ONE, ZERO + PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) +* .. +* .. Local Scalars .. + LOGICAL UPPER + INTEGER J, KLD, KN + DOUBLE PRECISION AJJ +* .. +* .. External Functions .. + LOGICAL LSAME + EXTERNAL LSAME +* .. +* .. External Subroutines .. + EXTERNAL DSCAL, DSYR, XERBLA +* .. +* .. Intrinsic Functions .. + INTRINSIC MAX, MIN, SQRT +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* + INFO = 0 + UPPER = LSAME( UPLO, 'U' ) + IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN + INFO = -1 + ELSE IF( N.LT.0 ) THEN + INFO = -2 + ELSE IF( KD.LT.0 ) THEN + INFO = -3 + ELSE IF( LDAB.LT.KD+1 ) THEN + INFO = -5 + END IF + IF( INFO.NE.0 ) THEN + CALL XERBLA( 'DPBTF2', -INFO ) + RETURN + END IF +* +* Quick return if possible +* + IF( N.EQ.0 ) + $ RETURN +* + KLD = MAX( 1, LDAB-1 ) +* + IF( UPPER ) THEN +* +* Compute the Cholesky factorization A = U'*U. +* + DO 10 J = 1, N +* +* Compute U(J,J) and test for non-positive-definiteness. +* + AJJ = AB( KD+1, J ) + IF( AJJ.LE.ZERO ) + $ GO TO 30 + AJJ = SQRT( AJJ ) + AB( KD+1, J ) = AJJ +* +* Compute elements J+1:J+KN of row J and update the +* trailing submatrix within the band. +* + KN = MIN( KD, N-J ) + IF( KN.GT.0 ) THEN + CALL DSCAL( KN, ONE / AJJ, AB( KD, J+1 ), KLD ) + CALL DSYR( 'Upper', KN, -ONE, AB( KD, J+1 ), KLD, + $ AB( KD+1, J+1 ), KLD ) + END IF + 10 CONTINUE + ELSE +* +* Compute the Cholesky factorization A = L*L'. +* + DO 20 J = 1, N +* +* Compute L(J,J) and test for non-positive-definiteness. +* + AJJ = AB( 1, J ) + IF( AJJ.LE.ZERO ) + $ GO TO 30 + AJJ = SQRT( AJJ ) + AB( 1, J ) = AJJ +* +* Compute elements J+1:J+KN of column J and update the +* trailing submatrix within the band. +* + KN = MIN( KD, N-J ) + IF( KN.GT.0 ) THEN + CALL DSCAL( KN, ONE / AJJ, AB( 2, J ), 1 ) + CALL DSYR( 'Lower', KN, -ONE, AB( 2, J ), 1, + $ AB( 1, J+1 ), KLD ) + END IF + 20 CONTINUE + END IF + RETURN +* + 30 CONTINUE + INFO = J + RETURN +* +* End of DPBTF2 +* + END + LOGICAL FUNCTION LSAME( CA, CB ) +* +* -- LAPACK auxiliary routine (version 2.0) -- +* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., +* Courant Institute, Argonne National Lab, and Rice University +* September 30, 1994 +* +* .. Scalar Arguments .. + CHARACTER CA, CB +* .. +* +* Purpose +* ======= +* +* LSAME returns .TRUE. if CA is the same letter as CB regardless of +* case. +* +* Arguments +* ========= +* +* CA (input) CHARACTER*1 +* CB (input) CHARACTER*1 +* CA and CB specify the single characters to be compared. +* +* ===================================================================== +* +* .. Intrinsic Functions .. + INTRINSIC ICHAR +* .. +* .. Local Scalars .. + INTEGER INTA, INTB, ZCODE +* .. +* .. Executable Statements .. +* +* Test if the characters are equal +* + LSAME = CA.EQ.CB + IF( LSAME ) + $ RETURN +* +* Now test for equivalence if both characters are alphabetic. +* + ZCODE = ICHAR( 'Z' ) +* +* Use 'Z' rather than 'A' so that ASCII can be detected on Prime +* machines, on which ICHAR returns a value with bit 8 set. +* ICHAR('A') on Prime machines returns 193 which is the same as +* ICHAR('A') on an EBCDIC machine. +* + INTA = ICHAR( CA ) + INTB = ICHAR( CB ) +* + IF( ZCODE.EQ.90 .OR. ZCODE.EQ.122 ) THEN +* +* ASCII is assumed - ZCODE is the ASCII code of either lower or +* upper case 'Z'. +* + IF( INTA.GE.97 .AND. INTA.LE.122 ) INTA = INTA - 32 + IF( INTB.GE.97 .AND. INTB.LE.122 ) INTB = INTB - 32 +* + ELSE IF( ZCODE.EQ.233 .OR. ZCODE.EQ.169 ) THEN +* +* EBCDIC is assumed - ZCODE is the EBCDIC code of either lower or +* upper case 'Z'. +* + IF( INTA.GE.129 .AND. INTA.LE.137 .OR. + $ INTA.GE.145 .AND. INTA.LE.153 .OR. + $ INTA.GE.162 .AND. INTA.LE.169 ) INTA = INTA + 64 + IF( INTB.GE.129 .AND. INTB.LE.137 .OR. + $ INTB.GE.145 .AND. INTB.LE.153 .OR. + $ INTB.GE.162 .AND. INTB.LE.169 ) INTB = INTB + 64 +* + ELSE IF( ZCODE.EQ.218 .OR. ZCODE.EQ.250 ) THEN +* +* ASCII is assumed, on Prime machines - ZCODE is the ASCII code +* plus 128 of either lower or upper case 'Z'. +* + IF( INTA.GE.225 .AND. INTA.LE.250 ) INTA = INTA - 32 + IF( INTB.GE.225 .AND. INTB.LE.250 ) INTB = INTB - 32 + END IF + LSAME = INTA.EQ.INTB +* +* RETURN +* +* End of LSAME +* + END + SUBROUTINE XERBLA( SRNAME, INFO ) +* +* -- LAPACK auxiliary routine (version 2.0) -- +* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., +* Courant Institute, Argonne National Lab, and Rice University +* September 30, 1994 +* +* .. Scalar Arguments .. + CHARACTER*6 SRNAME + INTEGER INFO +* .. +* +* Purpose +* ======= +* +* XERBLA is an error handler for the LAPACK routines. +* It is called by an LAPACK routine if an input parameter has an +* invalid value. A message is printed and execution stops. +* +* Installers may consider modifying the STOP statement in order to +* call system-specific exception-handling facilities. +* +* Arguments +* ========= +* +* SRNAME (input) CHARACTER*6 +* The name of the routine which called XERBLA. +* +* INFO (input) INTEGER +* The position of the invalid parameter in the parameter list +* of the calling routine. +* +* ===================================================================== +* +* .. Executable Statements .. +* + WRITE( *, FMT = 9999 )SRNAME, INFO +* + STOP +* + 9999 FORMAT( ' ** On entry to ', A6, ' parameter number ', I2, ' had ', + $ 'an illegal value' ) +* +* End of XERBLA +* + END + SUBROUTINE DPOTF2( UPLO, N, A, LDA, INFO ) +* +* -- LAPACK routine (version 2.0) -- +* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., +* Courant Institute, Argonne National Lab, and Rice University +* February 29, 1992 +* +* .. Scalar Arguments .. + CHARACTER UPLO + INTEGER INFO, LDA, N +* .. +* .. Array Arguments .. + DOUBLE PRECISION A( LDA, * ) +* .. +* +* Purpose +* ======= +* +* DPOTF2 computes the Cholesky factorization of a real symmetric +* positive definite matrix A. +* +* The factorization has the form +* A = U' * U , if UPLO = 'U', or +* A = L * L', if UPLO = 'L', +* where U is an upper triangular matrix and L is lower triangular. +* +* This is the unblocked version of the algorithm, calling Level 2 BLAS. +* +* Arguments +* ========= +* +* UPLO (input) CHARACTER*1 +* Specifies whether the upper or lower triangular part of the +* symmetric matrix A is stored. +* = 'U': Upper triangular +* = 'L': Lower triangular +* +* N (input) INTEGER +* The order of the matrix A. N >= 0. +* +* A (input/output) DOUBLE PRECISION array, dimension (LDA,N) +* On entry, the symmetric matrix A. If UPLO = 'U', the leading +* n by n upper triangular part of A contains the upper +* triangular part of the matrix A, and the strictly lower +* triangular part of A is not referenced. If UPLO = 'L', the +* leading n by n lower triangular part of A contains the lower +* triangular part of the matrix A, and the strictly upper +* triangular part of A is not referenced. +* +* On exit, if INFO = 0, the factor U or L from the Cholesky +* factorization A = U'*U or A = L*L'. +* +* LDA (input) INTEGER +* The leading dimension of the array A. LDA >= max(1,N). +* +* INFO (output) INTEGER +* = 0: successful exit +* < 0: if INFO = -k, the k-th argument had an illegal value +* > 0: if INFO = k, the leading minor of order k is not +* positive definite, and the factorization could not be +* completed. +* +* ===================================================================== +* +* .. Parameters .. + DOUBLE PRECISION ONE, ZERO + PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) +* .. +* .. Local Scalars .. + LOGICAL UPPER + INTEGER J + DOUBLE PRECISION AJJ +* .. +* .. External Functions .. + LOGICAL LSAME + DOUBLE PRECISION DDOT + EXTERNAL LSAME, DDOT +* .. +* .. External Subroutines .. + EXTERNAL DGEMV, DSCAL, XERBLA +* .. +* .. Intrinsic Functions .. + INTRINSIC MAX, SQRT +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* + INFO = 0 + UPPER = LSAME( UPLO, 'U' ) + IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN + INFO = -1 + ELSE IF( N.LT.0 ) THEN + INFO = -2 + ELSE IF( LDA.LT.MAX( 1, N ) ) THEN + INFO = -4 + END IF + IF( INFO.NE.0 ) THEN + CALL XERBLA( 'DPOTF2', -INFO ) + RETURN + END IF +* +* Quick return if possible +* + IF( N.EQ.0 ) + $ RETURN +* + IF( UPPER ) THEN +* +* Compute the Cholesky factorization A = U'*U. +* + DO 10 J = 1, N +* +* Compute U(J,J) and test for non-positive-definiteness. +* + AJJ = A( J, J ) - DDOT( J-1, A( 1, J ), 1, A( 1, J ), 1 ) + IF( AJJ.LE.ZERO ) THEN + A( J, J ) = AJJ + GO TO 30 + END IF + AJJ = SQRT( AJJ ) + A( J, J ) = AJJ +* +* Compute elements J+1:N of row J. +* + IF( J.LT.N ) THEN + CALL DGEMV( 'Transpose', J-1, N-J, -ONE, A( 1, J+1 ), + $ LDA, A( 1, J ), 1, ONE, A( J, J+1 ), LDA ) + CALL DSCAL( N-J, ONE / AJJ, A( J, J+1 ), LDA ) + END IF + 10 CONTINUE + ELSE +* +* Compute the Cholesky factorization A = L*L'. +* + DO 20 J = 1, N +* +* Compute L(J,J) and test for non-positive-definiteness. +* + AJJ = A( J, J ) - DDOT( J-1, A( J, 1 ), LDA, A( J, 1 ), + $ LDA ) + IF( AJJ.LE.ZERO ) THEN + A( J, J ) = AJJ + GO TO 30 + END IF + AJJ = SQRT( AJJ ) + A( J, J ) = AJJ +* +* Compute elements J+1:N of column J. +* + IF( J.LT.N ) THEN + CALL DGEMV( 'No transpose', N-J, J-1, -ONE, A( J+1, 1 ), + $ LDA, A( J, 1 ), LDA, ONE, A( J+1, J ), 1 ) + CALL DSCAL( N-J, ONE / AJJ, A( J+1, J ), 1 ) + END IF + 20 CONTINUE + END IF + GO TO 40 +* + 30 CONTINUE + INFO = J +* + 40 CONTINUE + RETURN +* +* End of DPOTF2 +* + END + subroutine dscal(n,da,dx,incx) +c +c scales a vector by a constant. +c uses unrolled loops for increment equal to one. +c jack dongarra, linpack, 3/11/78. +c modified 3/93 to return if incx .le. 0. +c modified 12/3/93, array(1) declarations changed to array(*) +c + double precision da,dx(*) + integer i,incx,m,mp1,n,nincx +c + if( n.le.0 .or. incx.le.0 )return + if(incx.eq.1)go to 20 +c +c code for increment not equal to 1 +c + nincx = n*incx + do 10 i = 1,nincx,incx + dx(i) = da*dx(i) + 10 continue + return +c +c code for increment equal to 1 +c +c +c clean-up loop +c + 20 m = mod(n,5) + if( m .eq. 0 ) go to 40 + do 30 i = 1,m + dx(i) = da*dx(i) + 30 continue + if( n .lt. 5 ) return + 40 mp1 = m + 1 + do 50 i = mp1,n,5 + dx(i) = da*dx(i) + dx(i + 1) = da*dx(i + 1) + dx(i + 2) = da*dx(i + 2) + dx(i + 3) = da*dx(i + 3) + dx(i + 4) = da*dx(i + 4) + 50 continue + return + end + double precision function ddot(n,dx,incx,dy,incy) +c +c forms the dot product of two vectors. +c uses unrolled loops for increments equal to one. +c jack dongarra, linpack, 3/11/78. +c modified 12/3/93, array(1) declarations changed to array(*) +c + double precision dx(*),dy(*),dtemp + integer i,incx,incy,ix,iy,m,mp1,n +c + ddot = 0.0d0 + dtemp = 0.0d0 + if(n.le.0)return + if(incx.eq.1.and.incy.eq.1)go to 20 +c +c code for unequal increments or equal increments +c not equal to 1 +c + ix = 1 + iy = 1 + if(incx.lt.0)ix = (-n+1)*incx + 1 + if(incy.lt.0)iy = (-n+1)*incy + 1 + do 10 i = 1,n + dtemp = dtemp + dx(ix)*dy(iy) + ix = ix + incx + iy = iy + incy + 10 continue + ddot = dtemp + return +c +c code for both increments equal to 1 +c +c +c clean-up loop +c + 20 m = mod(n,5) + if( m .eq. 0 ) go to 40 + do 30 i = 1,m + dtemp = dtemp + dx(i)*dy(i) + 30 continue + if( n .lt. 5 ) go to 60 + 40 mp1 = m + 1 + do 50 i = mp1,n,5 + dtemp = dtemp + dx(i)*dy(i) + dx(i + 1)*dy(i + 1) + + * dx(i + 2)*dy(i + 2) + dx(i + 3)*dy(i + 3) + dx(i + 4)*dy(i + 4) + 50 continue + 60 ddot = dtemp + return + end + SUBROUTINE DTRSM ( SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, + $ B, LDB ) +* .. Scalar Arguments .. + CHARACTER*1 SIDE, UPLO, TRANSA, DIAG + INTEGER M, N, LDA, LDB + DOUBLE PRECISION ALPHA +* .. Array Arguments .. + DOUBLE PRECISION A( LDA, * ), B( LDB, * ) +* .. +* +* Purpose +* ======= +* +* DTRSM solves one of the matrix equations +* +* op( A )*X = alpha*B, or X*op( A ) = alpha*B, +* +* where alpha is a scalar, X and B are m by n matrices, A is a unit, or +* non-unit, upper or lower triangular matrix and op( A ) is one of +* +* op( A ) = A or op( A ) = A'. +* +* The matrix X is overwritten on B. +* +* Parameters +* ========== +* +* SIDE - CHARACTER*1. +* On entry, SIDE specifies whether op( A ) appears on the left +* or right of X as follows: +* +* SIDE = 'L' or 'l' op( A )*X = alpha*B. +* +* SIDE = 'R' or 'r' X*op( A ) = alpha*B. +* +* Unchanged on exit. +* +* UPLO - CHARACTER*1. +* On entry, UPLO specifies whether the matrix A is an upper or +* lower triangular matrix as follows: +* +* UPLO = 'U' or 'u' A is an upper triangular matrix. +* +* UPLO = 'L' or 'l' A is a lower triangular matrix. +* +* Unchanged on exit. +* +* TRANSA - CHARACTER*1. +* On entry, TRANSA specifies the form of op( A ) to be used in +* the matrix multiplication as follows: +* +* TRANSA = 'N' or 'n' op( A ) = A. +* +* TRANSA = 'T' or 't' op( A ) = A'. +* +* TRANSA = 'C' or 'c' op( A ) = A'. +* +* Unchanged on exit. +* +* DIAG - CHARACTER*1. +* On entry, DIAG specifies whether or not A is unit triangular +* as follows: +* +* DIAG = 'U' or 'u' A is assumed to be unit triangular. +* +* DIAG = 'N' or 'n' A is not assumed to be unit +* triangular. +* +* Unchanged on exit. +* +* M - INTEGER. +* On entry, M specifies the number of rows of B. M must be at +* least zero. +* Unchanged on exit. +* +* N - INTEGER. +* On entry, N specifies the number of columns of B. N must be +* at least zero. +* Unchanged on exit. +* +* ALPHA - DOUBLE PRECISION. +* On entry, ALPHA specifies the scalar alpha. When alpha is +* zero then A is not referenced and B need not be set before +* entry. +* Unchanged on exit. +* +* A - DOUBLE PRECISION array of DIMENSION ( LDA, k ), where k is m +* when SIDE = 'L' or 'l' and is n when SIDE = 'R' or 'r'. +* Before entry with UPLO = 'U' or 'u', the leading k by k +* upper triangular part of the array A must contain the upper +* triangular matrix and the strictly lower triangular part of +* A is not referenced. +* Before entry with UPLO = 'L' or 'l', the leading k by k +* lower triangular part of the array A must contain the lower +* triangular matrix and the strictly upper triangular part of +* A is not referenced. +* Note that when DIAG = 'U' or 'u', the diagonal elements of +* A are not referenced either, but are assumed to be unity. +* Unchanged on exit. +* +* LDA - INTEGER. +* On entry, LDA specifies the first dimension of A as declared +* in the calling (sub) program. When SIDE = 'L' or 'l' then +* LDA must be at least max( 1, m ), when SIDE = 'R' or 'r' +* then LDA must be at least max( 1, n ). +* Unchanged on exit. +* +* B - DOUBLE PRECISION array of DIMENSION ( LDB, n ). +* Before entry, the leading m by n part of the array B must +* contain the right-hand side matrix B, and on exit is +* overwritten by the solution matrix X. +* +* LDB - INTEGER. +* On entry, LDB specifies the first dimension of B as declared +* in the calling (sub) program. LDB must be at least +* max( 1, m ). +* Unchanged on exit. +* +* +* Level 3 Blas routine. +* +* +* -- Written on 8-February-1989. +* Jack Dongarra, Argonne National Laboratory. +* Iain Duff, AERE Harwell. +* Jeremy Du Croz, Numerical Algorithms Group Ltd. +* Sven Hammarling, Numerical Algorithms Group Ltd. +* +* +* .. External Functions .. + LOGICAL LSAME + EXTERNAL LSAME +* .. External Subroutines .. + EXTERNAL XERBLA +* .. Intrinsic Functions .. + INTRINSIC MAX +* .. Local Scalars .. + LOGICAL LSIDE, NOUNIT, UPPER + INTEGER I, INFO, J, K, NROWA + DOUBLE PRECISION TEMP +* .. Parameters .. + DOUBLE PRECISION ONE , ZERO + PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* + LSIDE = LSAME( SIDE , 'L' ) + IF( LSIDE )THEN + NROWA = M + ELSE + NROWA = N + END IF + NOUNIT = LSAME( DIAG , 'N' ) + UPPER = LSAME( UPLO , 'U' ) +* + INFO = 0 + IF( ( .NOT.LSIDE ).AND. + $ ( .NOT.LSAME( SIDE , 'R' ) ) )THEN + INFO = 1 + ELSE IF( ( .NOT.UPPER ).AND. + $ ( .NOT.LSAME( UPLO , 'L' ) ) )THEN + INFO = 2 + ELSE IF( ( .NOT.LSAME( TRANSA, 'N' ) ).AND. + $ ( .NOT.LSAME( TRANSA, 'T' ) ).AND. + $ ( .NOT.LSAME( TRANSA, 'C' ) ) )THEN + INFO = 3 + ELSE IF( ( .NOT.LSAME( DIAG , 'U' ) ).AND. + $ ( .NOT.LSAME( DIAG , 'N' ) ) )THEN + INFO = 4 + ELSE IF( M .LT.0 )THEN + INFO = 5 + ELSE IF( N .LT.0 )THEN + INFO = 6 + ELSE IF( LDA.LT.MAX( 1, NROWA ) )THEN + INFO = 9 + ELSE IF( LDB.LT.MAX( 1, M ) )THEN + INFO = 11 + END IF + IF( INFO.NE.0 )THEN + CALL XERBLA( 'DTRSM ', INFO ) + RETURN + END IF +* +* Quick return if possible. +* + IF( N.EQ.0 ) + $ RETURN +* +* And when alpha.eq.zero. +* + IF( ALPHA.EQ.ZERO )THEN + DO 20, J = 1, N + DO 10, I = 1, M + B( I, J ) = ZERO + 10 CONTINUE + 20 CONTINUE + RETURN + END IF +* +* Start the operations. +* + IF( LSIDE )THEN + IF( LSAME( TRANSA, 'N' ) )THEN +* +* Form B := alpha*inv( A )*B. +* + IF( UPPER )THEN + DO 60, J = 1, N + IF( ALPHA.NE.ONE )THEN + DO 30, I = 1, M + B( I, J ) = ALPHA*B( I, J ) + 30 CONTINUE + END IF + DO 50, K = M, 1, -1 + IF( B( K, J ).NE.ZERO )THEN + IF( NOUNIT ) + $ B( K, J ) = B( K, J )/A( K, K ) + DO 40, I = 1, K - 1 + B( I, J ) = B( I, J ) - B( K, J )*A( I, K ) + 40 CONTINUE + END IF + 50 CONTINUE + 60 CONTINUE + ELSE + DO 100, J = 1, N + IF( ALPHA.NE.ONE )THEN + DO 70, I = 1, M + B( I, J ) = ALPHA*B( I, J ) + 70 CONTINUE + END IF + DO 90 K = 1, M + IF( B( K, J ).NE.ZERO )THEN + IF( NOUNIT ) + $ B( K, J ) = B( K, J )/A( K, K ) + DO 80, I = K + 1, M + B( I, J ) = B( I, J ) - B( K, J )*A( I, K ) + 80 CONTINUE + END IF + 90 CONTINUE + 100 CONTINUE + END IF + ELSE +* +* Form B := alpha*inv( A' )*B. +* + IF( UPPER )THEN + DO 130, J = 1, N + DO 120, I = 1, M + TEMP = ALPHA*B( I, J ) + DO 110, K = 1, I - 1 + TEMP = TEMP - A( K, I )*B( K, J ) + 110 CONTINUE + IF( NOUNIT ) + $ TEMP = TEMP/A( I, I ) + B( I, J ) = TEMP + 120 CONTINUE + 130 CONTINUE + ELSE + DO 160, J = 1, N + DO 150, I = M, 1, -1 + TEMP = ALPHA*B( I, J ) + DO 140, K = I + 1, M + TEMP = TEMP - A( K, I )*B( K, J ) + 140 CONTINUE + IF( NOUNIT ) + $ TEMP = TEMP/A( I, I ) + B( I, J ) = TEMP + 150 CONTINUE + 160 CONTINUE + END IF + END IF + ELSE + IF( LSAME( TRANSA, 'N' ) )THEN +* +* Form B := alpha*B*inv( A ). +* + IF( UPPER )THEN + DO 210, J = 1, N + IF( ALPHA.NE.ONE )THEN + DO 170, I = 1, M + B( I, J ) = ALPHA*B( I, J ) + 170 CONTINUE + END IF + DO 190, K = 1, J - 1 + IF( A( K, J ).NE.ZERO )THEN + DO 180, I = 1, M + B( I, J ) = B( I, J ) - A( K, J )*B( I, K ) + 180 CONTINUE + END IF + 190 CONTINUE + IF( NOUNIT )THEN + TEMP = ONE/A( J, J ) + DO 200, I = 1, M + B( I, J ) = TEMP*B( I, J ) + 200 CONTINUE + END IF + 210 CONTINUE + ELSE + DO 260, J = N, 1, -1 + IF( ALPHA.NE.ONE )THEN + DO 220, I = 1, M + B( I, J ) = ALPHA*B( I, J ) + 220 CONTINUE + END IF + DO 240, K = J + 1, N + IF( A( K, J ).NE.ZERO )THEN + DO 230, I = 1, M + B( I, J ) = B( I, J ) - A( K, J )*B( I, K ) + 230 CONTINUE + END IF + 240 CONTINUE + IF( NOUNIT )THEN + TEMP = ONE/A( J, J ) + DO 250, I = 1, M + B( I, J ) = TEMP*B( I, J ) + 250 CONTINUE + END IF + 260 CONTINUE + END IF + ELSE +* +* Form B := alpha*B*inv( A' ). +* + IF( UPPER )THEN + DO 310, K = N, 1, -1 + IF( NOUNIT )THEN + TEMP = ONE/A( K, K ) + DO 270, I = 1, M + B( I, K ) = TEMP*B( I, K ) + 270 CONTINUE + END IF + DO 290, J = 1, K - 1 + IF( A( J, K ).NE.ZERO )THEN + TEMP = A( J, K ) + DO 280, I = 1, M + B( I, J ) = B( I, J ) - TEMP*B( I, K ) + 280 CONTINUE + END IF + 290 CONTINUE + IF( ALPHA.NE.ONE )THEN + DO 300, I = 1, M + B( I, K ) = ALPHA*B( I, K ) + 300 CONTINUE + END IF + 310 CONTINUE + ELSE + DO 360, K = 1, N + IF( NOUNIT )THEN + TEMP = ONE/A( K, K ) + DO 320, I = 1, M + B( I, K ) = TEMP*B( I, K ) + 320 CONTINUE + END IF + DO 340, J = K + 1, N + IF( A( J, K ).NE.ZERO )THEN + TEMP = A( J, K ) + DO 330, I = 1, M + B( I, J ) = B( I, J ) - TEMP*B( I, K ) + 330 CONTINUE + END IF + 340 CONTINUE + IF( ALPHA.NE.ONE )THEN + DO 350, I = 1, M + B( I, K ) = ALPHA*B( I, K ) + 350 CONTINUE + END IF + 360 CONTINUE + END IF + END IF + END IF +* + RETURN +* +* End of DTRSM . +* + END + SUBROUTINE DTBSV ( UPLO, TRANS, DIAG, N, K, A, LDA, X, INCX ) +* .. Scalar Arguments .. + INTEGER INCX, K, LDA, N + CHARACTER*1 DIAG, TRANS, UPLO +* .. Array Arguments .. + DOUBLE PRECISION A( LDA, * ), X( * ) +* .. +* +* Purpose +* ======= +* +* DTBSV solves one of the systems of equations +* +* A*x = b, or A'*x = b, +* +* where b and x are n element vectors and A is an n by n unit, or +* non-unit, upper or lower triangular band matrix, with ( k + 1 ) +* diagonals. +* +* No test for singularity or near-singularity is included in this +* routine. Such tests must be performed before calling this routine. +* +* Parameters +* ========== +* +* UPLO - CHARACTER*1. +* On entry, UPLO specifies whether the matrix is an upper or +* lower triangular matrix as follows: +* +* UPLO = 'U' or 'u' A is an upper triangular matrix. +* +* UPLO = 'L' or 'l' A is a lower triangular matrix. +* +* Unchanged on exit. +* +* TRANS - CHARACTER*1. +* On entry, TRANS specifies the equations to be solved as +* follows: +* +* TRANS = 'N' or 'n' A*x = b. +* +* TRANS = 'T' or 't' A'*x = b. +* +* TRANS = 'C' or 'c' A'*x = b. +* +* Unchanged on exit. +* +* DIAG - CHARACTER*1. +* On entry, DIAG specifies whether or not A is unit +* triangular as follows: +* +* DIAG = 'U' or 'u' A is assumed to be unit triangular. +* +* DIAG = 'N' or 'n' A is not assumed to be unit +* triangular. +* +* Unchanged on exit. +* +* N - INTEGER. +* On entry, N specifies the order of the matrix A. +* N must be at least zero. +* Unchanged on exit. +* +* K - INTEGER. +* On entry with UPLO = 'U' or 'u', K specifies the number of +* super-diagonals of the matrix A. +* On entry with UPLO = 'L' or 'l', K specifies the number of +* sub-diagonals of the matrix A. +* K must satisfy 0 .le. K. +* Unchanged on exit. +* +* A - DOUBLE PRECISION array of DIMENSION ( LDA, n ). +* Before entry with UPLO = 'U' or 'u', the leading ( k + 1 ) +* by n part of the array A must contain the upper triangular +* band part of the matrix of coefficients, supplied column by +* column, with the leading diagonal of the matrix in row +* ( k + 1 ) of the array, the first super-diagonal starting at +* position 2 in row k, and so on. The top left k by k triangle +* of the array A is not referenced. +* The following program segment will transfer an upper +* triangular band matrix from conventional full matrix storage +* to band storage: +* +* DO 20, J = 1, N +* M = K + 1 - J +* DO 10, I = MAX( 1, J - K ), J +* A( M + I, J ) = matrix( I, J ) +* 10 CONTINUE +* 20 CONTINUE +* +* Before entry with UPLO = 'L' or 'l', the leading ( k + 1 ) +* by n part of the array A must contain the lower triangular +* band part of the matrix of coefficients, supplied column by +* column, with the leading diagonal of the matrix in row 1 of +* the array, the first sub-diagonal starting at position 1 in +* row 2, and so on. The bottom right k by k triangle of the +* array A is not referenced. +* The following program segment will transfer a lower +* triangular band matrix from conventional full matrix storage +* to band storage: +* +* DO 20, J = 1, N +* M = 1 - J +* DO 10, I = J, MIN( N, J + K ) +* A( M + I, J ) = matrix( I, J ) +* 10 CONTINUE +* 20 CONTINUE +* +* Note that when DIAG = 'U' or 'u' the elements of the array A +* corresponding to the diagonal elements of the matrix are not +* referenced, but are assumed to be unity. +* Unchanged on exit. +* +* LDA - INTEGER. +* On entry, LDA specifies the first dimension of A as declared +* in the calling (sub) program. LDA must be at least +* ( k + 1 ). +* Unchanged on exit. +* +* X - DOUBLE PRECISION array of dimension at least +* ( 1 + ( n - 1 )*abs( INCX ) ). +* Before entry, the incremented array X must contain the n +* element right-hand side vector b. On exit, X is overwritten +* with the solution vector x. +* +* INCX - INTEGER. +* On entry, INCX specifies the increment for the elements of +* X. INCX must not be zero. +* Unchanged on exit. +* +* +* Level 2 Blas routine. +* +* -- Written on 22-October-1986. +* Jack Dongarra, Argonne National Lab. +* Jeremy Du Croz, Nag Central Office. +* Sven Hammarling, Nag Central Office. +* Richard Hanson, Sandia National Labs. +* +* +* .. Parameters .. + DOUBLE PRECISION ZERO + PARAMETER ( ZERO = 0.0D+0 ) +* .. Local Scalars .. + DOUBLE PRECISION TEMP + INTEGER I, INFO, IX, J, JX, KPLUS1, KX, L + LOGICAL NOUNIT +* .. External Functions .. + LOGICAL LSAME + EXTERNAL LSAME +* .. External Subroutines .. + EXTERNAL XERBLA +* .. Intrinsic Functions .. + INTRINSIC MAX, MIN +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* + INFO = 0 + IF ( .NOT.LSAME( UPLO , 'U' ).AND. + $ .NOT.LSAME( UPLO , 'L' ) )THEN + INFO = 1 + ELSE IF( .NOT.LSAME( TRANS, 'N' ).AND. + $ .NOT.LSAME( TRANS, 'T' ).AND. + $ .NOT.LSAME( TRANS, 'C' ) )THEN + INFO = 2 + ELSE IF( .NOT.LSAME( DIAG , 'U' ).AND. + $ .NOT.LSAME( DIAG , 'N' ) )THEN + INFO = 3 + ELSE IF( N.LT.0 )THEN + INFO = 4 + ELSE IF( K.LT.0 )THEN + INFO = 5 + ELSE IF( LDA.LT.( K + 1 ) )THEN + INFO = 7 + ELSE IF( INCX.EQ.0 )THEN + INFO = 9 + END IF + IF( INFO.NE.0 )THEN + CALL XERBLA( 'DTBSV ', INFO ) + RETURN + END IF +* +* Quick return if possible. +* + IF( N.EQ.0 ) + $ RETURN +* + NOUNIT = LSAME( DIAG, 'N' ) +* +* Set up the start point in X if the increment is not unity. This +* will be ( N - 1 )*INCX too small for descending loops. +* + IF( INCX.LE.0 )THEN + KX = 1 - ( N - 1 )*INCX + ELSE IF( INCX.NE.1 )THEN + KX = 1 + END IF +* +* Start the operations. In this version the elements of A are +* accessed by sequentially with one pass through A. +* + IF( LSAME( TRANS, 'N' ) )THEN +* +* Form x := inv( A )*x. +* + IF( LSAME( UPLO, 'U' ) )THEN + KPLUS1 = K + 1 + IF( INCX.EQ.1 )THEN + DO 20, J = N, 1, -1 + IF( X( J ).NE.ZERO )THEN + L = KPLUS1 - J + IF( NOUNIT ) + $ X( J ) = X( J )/A( KPLUS1, J ) + TEMP = X( J ) + DO 10, I = J - 1, MAX( 1, J - K ), -1 + X( I ) = X( I ) - TEMP*A( L + I, J ) + 10 CONTINUE + END IF + 20 CONTINUE + ELSE + KX = KX + ( N - 1 )*INCX + JX = KX + DO 40, J = N, 1, -1 + KX = KX - INCX + IF( X( JX ).NE.ZERO )THEN + IX = KX + L = KPLUS1 - J + IF( NOUNIT ) + $ X( JX ) = X( JX )/A( KPLUS1, J ) + TEMP = X( JX ) + DO 30, I = J - 1, MAX( 1, J - K ), -1 + X( IX ) = X( IX ) - TEMP*A( L + I, J ) + IX = IX - INCX + 30 CONTINUE + END IF + JX = JX - INCX + 40 CONTINUE + END IF + ELSE + IF( INCX.EQ.1 )THEN + DO 60, J = 1, N + IF( X( J ).NE.ZERO )THEN + L = 1 - J + IF( NOUNIT ) + $ X( J ) = X( J )/A( 1, J ) + TEMP = X( J ) + DO 50, I = J + 1, MIN( N, J + K ) + X( I ) = X( I ) - TEMP*A( L + I, J ) + 50 CONTINUE + END IF + 60 CONTINUE + ELSE + JX = KX + DO 80, J = 1, N + KX = KX + INCX + IF( X( JX ).NE.ZERO )THEN + IX = KX + L = 1 - J + IF( NOUNIT ) + $ X( JX ) = X( JX )/A( 1, J ) + TEMP = X( JX ) + DO 70, I = J + 1, MIN( N, J + K ) + X( IX ) = X( IX ) - TEMP*A( L + I, J ) + IX = IX + INCX + 70 CONTINUE + END IF + JX = JX + INCX + 80 CONTINUE + END IF + END IF + ELSE +* +* Form x := inv( A')*x. +* + IF( LSAME( UPLO, 'U' ) )THEN + KPLUS1 = K + 1 + IF( INCX.EQ.1 )THEN + DO 100, J = 1, N + TEMP = X( J ) + L = KPLUS1 - J + DO 90, I = MAX( 1, J - K ), J - 1 + TEMP = TEMP - A( L + I, J )*X( I ) + 90 CONTINUE + IF( NOUNIT ) + $ TEMP = TEMP/A( KPLUS1, J ) + X( J ) = TEMP + 100 CONTINUE + ELSE + JX = KX + DO 120, J = 1, N + TEMP = X( JX ) + IX = KX + L = KPLUS1 - J + DO 110, I = MAX( 1, J - K ), J - 1 + TEMP = TEMP - A( L + I, J )*X( IX ) + IX = IX + INCX + 110 CONTINUE + IF( NOUNIT ) + $ TEMP = TEMP/A( KPLUS1, J ) + X( JX ) = TEMP + JX = JX + INCX + IF( J.GT.K ) + $ KX = KX + INCX + 120 CONTINUE + END IF + ELSE + IF( INCX.EQ.1 )THEN + DO 140, J = N, 1, -1 + TEMP = X( J ) + L = 1 - J + DO 130, I = MIN( N, J + K ), J + 1, -1 + TEMP = TEMP - A( L + I, J )*X( I ) + 130 CONTINUE + IF( NOUNIT ) + $ TEMP = TEMP/A( 1, J ) + X( J ) = TEMP + 140 CONTINUE + ELSE + KX = KX + ( N - 1 )*INCX + JX = KX + DO 160, J = N, 1, -1 + TEMP = X( JX ) + IX = KX + L = 1 - J + DO 150, I = MIN( N, J + K ), J + 1, -1 + TEMP = TEMP - A( L + I, J )*X( IX ) + IX = IX - INCX + 150 CONTINUE + IF( NOUNIT ) + $ TEMP = TEMP/A( 1, J ) + X( JX ) = TEMP + JX = JX - INCX + IF( ( N - J ).GE.K ) + $ KX = KX - INCX + 160 CONTINUE + END IF + END IF + END IF +* + RETURN +* +* End of DTBSV . +* + END + SUBROUTINE DSYR ( UPLO, N, ALPHA, X, INCX, A, LDA ) +* .. Scalar Arguments .. + DOUBLE PRECISION ALPHA + INTEGER INCX, LDA, N + CHARACTER*1 UPLO +* .. Array Arguments .. + DOUBLE PRECISION A( LDA, * ), X( * ) +* .. +* +* Purpose +* ======= +* +* DSYR performs the symmetric rank 1 operation +* +* A := alpha*x*x' + A, +* +* where alpha is a real scalar, x is an n element vector and A is an +* n by n symmetric matrix. +* +* Parameters +* ========== +* +* UPLO - CHARACTER*1. +* On entry, UPLO specifies whether the upper or lower +* triangular part of the array A is to be referenced as +* follows: +* +* UPLO = 'U' or 'u' Only the upper triangular part of A +* is to be referenced. +* +* UPLO = 'L' or 'l' Only the lower triangular part of A +* is to be referenced. +* +* Unchanged on exit. +* +* N - INTEGER. +* On entry, N specifies the order of the matrix A. +* N must be at least zero. +* Unchanged on exit. +* +* ALPHA - DOUBLE PRECISION. +* On entry, ALPHA specifies the scalar alpha. +* Unchanged on exit. +* +* X - DOUBLE PRECISION array of dimension at least +* ( 1 + ( n - 1 )*abs( INCX ) ). +* Before entry, the incremented array X must contain the n +* element vector x. +* Unchanged on exit. +* +* INCX - INTEGER. +* On entry, INCX specifies the increment for the elements of +* X. INCX must not be zero. +* Unchanged on exit. +* +* A - DOUBLE PRECISION array of DIMENSION ( LDA, n ). +* Before entry with UPLO = 'U' or 'u', the leading n by n +* upper triangular part of the array A must contain the upper +* triangular part of the symmetric matrix and the strictly +* lower triangular part of A is not referenced. On exit, the +* upper triangular part of the array A is overwritten by the +* upper triangular part of the updated matrix. +* Before entry with UPLO = 'L' or 'l', the leading n by n +* lower triangular part of the array A must contain the lower +* triangular part of the symmetric matrix and the strictly +* upper triangular part of A is not referenced. On exit, the +* lower triangular part of the array A is overwritten by the +* lower triangular part of the updated matrix. +* +* LDA - INTEGER. +* On entry, LDA specifies the first dimension of A as declared +* in the calling (sub) program. LDA must be at least +* max( 1, n ). +* Unchanged on exit. +* +* +* Level 2 Blas routine. +* +* -- Written on 22-October-1986. +* Jack Dongarra, Argonne National Lab. +* Jeremy Du Croz, Nag Central Office. +* Sven Hammarling, Nag Central Office. +* Richard Hanson, Sandia National Labs. +* +* +* .. Parameters .. + DOUBLE PRECISION ZERO + PARAMETER ( ZERO = 0.0D+0 ) +* .. Local Scalars .. + DOUBLE PRECISION TEMP + INTEGER I, INFO, IX, J, JX, KX +* .. External Functions .. + LOGICAL LSAME + EXTERNAL LSAME +* .. External Subroutines .. + EXTERNAL XERBLA +* .. Intrinsic Functions .. + INTRINSIC MAX +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* + INFO = 0 + IF ( .NOT.LSAME( UPLO, 'U' ).AND. + $ .NOT.LSAME( UPLO, 'L' ) )THEN + INFO = 1 + ELSE IF( N.LT.0 )THEN + INFO = 2 + ELSE IF( INCX.EQ.0 )THEN + INFO = 5 + ELSE IF( LDA.LT.MAX( 1, N ) )THEN + INFO = 7 + END IF + IF( INFO.NE.0 )THEN + CALL XERBLA( 'DSYR ', INFO ) + RETURN + END IF +* +* Quick return if possible. +* + IF( ( N.EQ.0 ).OR.( ALPHA.EQ.ZERO ) ) + $ RETURN +* +* Set the start point in X if the increment is not unity. +* + IF( INCX.LE.0 )THEN + KX = 1 - ( N - 1 )*INCX + ELSE IF( INCX.NE.1 )THEN + KX = 1 + END IF +* +* Start the operations. In this version the elements of A are +* accessed sequentially with one pass through the triangular part +* of A. +* + IF( LSAME( UPLO, 'U' ) )THEN +* +* Form A when A is stored in upper triangle. +* + IF( INCX.EQ.1 )THEN + DO 20, J = 1, N + IF( X( J ).NE.ZERO )THEN + TEMP = ALPHA*X( J ) + DO 10, I = 1, J + A( I, J ) = A( I, J ) + X( I )*TEMP + 10 CONTINUE + END IF + 20 CONTINUE + ELSE + JX = KX + DO 40, J = 1, N + IF( X( JX ).NE.ZERO )THEN + TEMP = ALPHA*X( JX ) + IX = KX + DO 30, I = 1, J + A( I, J ) = A( I, J ) + X( IX )*TEMP + IX = IX + INCX + 30 CONTINUE + END IF + JX = JX + INCX + 40 CONTINUE + END IF + ELSE +* +* Form A when A is stored in lower triangle. +* + IF( INCX.EQ.1 )THEN + DO 60, J = 1, N + IF( X( J ).NE.ZERO )THEN + TEMP = ALPHA*X( J ) + DO 50, I = J, N + A( I, J ) = A( I, J ) + X( I )*TEMP + 50 CONTINUE + END IF + 60 CONTINUE + ELSE + JX = KX + DO 80, J = 1, N + IF( X( JX ).NE.ZERO )THEN + TEMP = ALPHA*X( JX ) + IX = JX + DO 70, I = J, N + A( I, J ) = A( I, J ) + X( IX )*TEMP + IX = IX + INCX + 70 CONTINUE + END IF + JX = JX + INCX + 80 CONTINUE + END IF + END IF +* + RETURN +* +* End of DSYR . +* + END + SUBROUTINE DSYRK ( UPLO, TRANS, N, K, ALPHA, A, LDA, + $ BETA, C, LDC ) +* .. Scalar Arguments .. + CHARACTER*1 UPLO, TRANS + INTEGER N, K, LDA, LDC + DOUBLE PRECISION ALPHA, BETA +* .. Array Arguments .. + DOUBLE PRECISION A( LDA, * ), C( LDC, * ) +* .. +* +* Purpose +* ======= +* +* DSYRK performs one of the symmetric rank k operations +* +* C := alpha*A*A' + beta*C, +* +* or +* +* C := alpha*A'*A + beta*C, +* +* where alpha and beta are scalars, C is an n by n symmetric matrix +* and A is an n by k matrix in the first case and a k by n matrix +* in the second case. +* +* Parameters +* ========== +* +* UPLO - CHARACTER*1. +* On entry, UPLO specifies whether the upper or lower +* triangular part of the array C is to be referenced as +* follows: +* +* UPLO = 'U' or 'u' Only the upper triangular part of C +* is to be referenced. +* +* UPLO = 'L' or 'l' Only the lower triangular part of C +* is to be referenced. +* +* Unchanged on exit. +* +* TRANS - CHARACTER*1. +* On entry, TRANS specifies the operation to be performed as +* follows: +* +* TRANS = 'N' or 'n' C := alpha*A*A' + beta*C. +* +* TRANS = 'T' or 't' C := alpha*A'*A + beta*C. +* +* TRANS = 'C' or 'c' C := alpha*A'*A + beta*C. +* +* Unchanged on exit. +* +* N - INTEGER. +* On entry, N specifies the order of the matrix C. N must be +* at least zero. +* Unchanged on exit. +* +* K - INTEGER. +* On entry with TRANS = 'N' or 'n', K specifies the number +* of columns of the matrix A, and on entry with +* TRANS = 'T' or 't' or 'C' or 'c', K specifies the number +* of rows of the matrix A. K must be at least zero. +* Unchanged on exit. +* +* ALPHA - DOUBLE PRECISION. +* On entry, ALPHA specifies the scalar alpha. +* Unchanged on exit. +* +* A - DOUBLE PRECISION array of DIMENSION ( LDA, ka ), where ka is +* k when TRANS = 'N' or 'n', and is n otherwise. +* Before entry with TRANS = 'N' or 'n', the leading n by k +* part of the array A must contain the matrix A, otherwise +* the leading k by n part of the array A must contain the +* matrix A. +* Unchanged on exit. +* +* LDA - INTEGER. +* On entry, LDA specifies the first dimension of A as declared +* in the calling (sub) program. When TRANS = 'N' or 'n' +* then LDA must be at least max( 1, n ), otherwise LDA must +* be at least max( 1, k ). +* Unchanged on exit. +* +* BETA - DOUBLE PRECISION. +* On entry, BETA specifies the scalar beta. +* Unchanged on exit. +* +* C - DOUBLE PRECISION array of DIMENSION ( LDC, n ). +* Before entry with UPLO = 'U' or 'u', the leading n by n +* upper triangular part of the array C must contain the upper +* triangular part of the symmetric matrix and the strictly +* lower triangular part of C is not referenced. On exit, the +* upper triangular part of the array C is overwritten by the +* upper triangular part of the updated matrix. +* Before entry with UPLO = 'L' or 'l', the leading n by n +* lower triangular part of the array C must contain the lower +* triangular part of the symmetric matrix and the strictly +* upper triangular part of C is not referenced. On exit, the +* lower triangular part of the array C is overwritten by the +* lower triangular part of the updated matrix. +* +* LDC - INTEGER. +* On entry, LDC specifies the first dimension of C as declared +* in the calling (sub) program. LDC must be at least +* max( 1, n ). +* Unchanged on exit. +* +* +* Level 3 Blas routine. +* +* -- Written on 8-February-1989. +* Jack Dongarra, Argonne National Laboratory. +* Iain Duff, AERE Harwell. +* Jeremy Du Croz, Numerical Algorithms Group Ltd. +* Sven Hammarling, Numerical Algorithms Group Ltd. +* +* +* .. External Functions .. + LOGICAL LSAME + EXTERNAL LSAME +* .. External Subroutines .. + EXTERNAL XERBLA +* .. Intrinsic Functions .. + INTRINSIC MAX +* .. Local Scalars .. + LOGICAL UPPER + INTEGER I, INFO, J, L, NROWA + DOUBLE PRECISION TEMP +* .. Parameters .. + DOUBLE PRECISION ONE , ZERO + PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* + IF( LSAME( TRANS, 'N' ) )THEN + NROWA = N + ELSE + NROWA = K + END IF + UPPER = LSAME( UPLO, 'U' ) +* + INFO = 0 + IF( ( .NOT.UPPER ).AND. + $ ( .NOT.LSAME( UPLO , 'L' ) ) )THEN + INFO = 1 + ELSE IF( ( .NOT.LSAME( TRANS, 'N' ) ).AND. + $ ( .NOT.LSAME( TRANS, 'T' ) ).AND. + $ ( .NOT.LSAME( TRANS, 'C' ) ) )THEN + INFO = 2 + ELSE IF( N .LT.0 )THEN + INFO = 3 + ELSE IF( K .LT.0 )THEN + INFO = 4 + ELSE IF( LDA.LT.MAX( 1, NROWA ) )THEN + INFO = 7 + ELSE IF( LDC.LT.MAX( 1, N ) )THEN + INFO = 10 + END IF + IF( INFO.NE.0 )THEN + CALL XERBLA( 'DSYRK ', INFO ) + RETURN + END IF +* +* Quick return if possible. +* + IF( ( N.EQ.0 ).OR. + $ ( ( ( ALPHA.EQ.ZERO ).OR.( K.EQ.0 ) ).AND.( BETA.EQ.ONE ) ) ) + $ RETURN +* +* And when alpha.eq.zero. +* + IF( ALPHA.EQ.ZERO )THEN + IF( UPPER )THEN + IF( BETA.EQ.ZERO )THEN + DO 20, J = 1, N + DO 10, I = 1, J + C( I, J ) = ZERO + 10 CONTINUE + 20 CONTINUE + ELSE + DO 40, J = 1, N + DO 30, I = 1, J + C( I, J ) = BETA*C( I, J ) + 30 CONTINUE + 40 CONTINUE + END IF + ELSE + IF( BETA.EQ.ZERO )THEN + DO 60, J = 1, N + DO 50, I = J, N + C( I, J ) = ZERO + 50 CONTINUE + 60 CONTINUE + ELSE + DO 80, J = 1, N + DO 70, I = J, N + C( I, J ) = BETA*C( I, J ) + 70 CONTINUE + 80 CONTINUE + END IF + END IF + RETURN + END IF +* +* Start the operations. +* + IF( LSAME( TRANS, 'N' ) )THEN +* +* Form C := alpha*A*A' + beta*C. +* + IF( UPPER )THEN + DO 130, J = 1, N + IF( BETA.EQ.ZERO )THEN + DO 90, I = 1, J + C( I, J ) = ZERO + 90 CONTINUE + ELSE IF( BETA.NE.ONE )THEN + DO 100, I = 1, J + C( I, J ) = BETA*C( I, J ) + 100 CONTINUE + END IF + DO 120, L = 1, K + IF( A( J, L ).NE.ZERO )THEN + TEMP = ALPHA*A( J, L ) + DO 110, I = 1, J + C( I, J ) = C( I, J ) + TEMP*A( I, L ) + 110 CONTINUE + END IF + 120 CONTINUE + 130 CONTINUE + ELSE + DO 180, J = 1, N + IF( BETA.EQ.ZERO )THEN + DO 140, I = J, N + C( I, J ) = ZERO + 140 CONTINUE + ELSE IF( BETA.NE.ONE )THEN + DO 150, I = J, N + C( I, J ) = BETA*C( I, J ) + 150 CONTINUE + END IF + DO 170, L = 1, K + IF( A( J, L ).NE.ZERO )THEN + TEMP = ALPHA*A( J, L ) + DO 160, I = J, N + C( I, J ) = C( I, J ) + TEMP*A( I, L ) + 160 CONTINUE + END IF + 170 CONTINUE + 180 CONTINUE + END IF + ELSE +* +* Form C := alpha*A'*A + beta*C. +* + IF( UPPER )THEN + DO 210, J = 1, N + DO 200, I = 1, J + TEMP = ZERO + DO 190, L = 1, K + TEMP = TEMP + A( L, I )*A( L, J ) + 190 CONTINUE + IF( BETA.EQ.ZERO )THEN + C( I, J ) = ALPHA*TEMP + ELSE + C( I, J ) = ALPHA*TEMP + BETA*C( I, J ) + END IF + 200 CONTINUE + 210 CONTINUE + ELSE + DO 240, J = 1, N + DO 230, I = J, N + TEMP = ZERO + DO 220, L = 1, K + TEMP = TEMP + A( L, I )*A( L, J ) + 220 CONTINUE + IF( BETA.EQ.ZERO )THEN + C( I, J ) = ALPHA*TEMP + ELSE + C( I, J ) = ALPHA*TEMP + BETA*C( I, J ) + END IF + 230 CONTINUE + 240 CONTINUE + END IF + END IF +* + RETURN +* +* End of DSYRK . +* + END + SUBROUTINE DGEMM ( TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, + $ BETA, C, LDC ) +* .. Scalar Arguments .. + CHARACTER*1 TRANSA, TRANSB + INTEGER M, N, K, LDA, LDB, LDC + DOUBLE PRECISION ALPHA, BETA +* .. Array Arguments .. + DOUBLE PRECISION A( LDA, * ), B( LDB, * ), C( LDC, * ) +* .. +* +* Purpose +* ======= +* +* DGEMM performs one of the matrix-matrix operations +* +* C := alpha*op( A )*op( B ) + beta*C, +* +* where op( X ) is one of +* +* op( X ) = X or op( X ) = X', +* +* alpha and beta are scalars, and A, B and C are matrices, with op( A ) +* an m by k matrix, op( B ) a k by n matrix and C an m by n matrix. +* +* Parameters +* ========== +* +* TRANSA - CHARACTER*1. +* On entry, TRANSA specifies the form of op( A ) to be used in +* the matrix multiplication as follows: +* +* TRANSA = 'N' or 'n', op( A ) = A. +* +* TRANSA = 'T' or 't', op( A ) = A'. +* +* TRANSA = 'C' or 'c', op( A ) = A'. +* +* Unchanged on exit. +* +* TRANSB - CHARACTER*1. +* On entry, TRANSB specifies the form of op( B ) to be used in +* the matrix multiplication as follows: +* +* TRANSB = 'N' or 'n', op( B ) = B. +* +* TRANSB = 'T' or 't', op( B ) = B'. +* +* TRANSB = 'C' or 'c', op( B ) = B'. +* +* Unchanged on exit. +* +* M - INTEGER. +* On entry, M specifies the number of rows of the matrix +* op( A ) and of the matrix C. M must be at least zero. +* Unchanged on exit. +* +* N - INTEGER. +* On entry, N specifies the number of columns of the matrix +* op( B ) and the number of columns of the matrix C. N must be +* at least zero. +* Unchanged on exit. +* +* K - INTEGER. +* On entry, K specifies the number of columns of the matrix +* op( A ) and the number of rows of the matrix op( B ). K must +* be at least zero. +* Unchanged on exit. +* +* ALPHA - DOUBLE PRECISION. +* On entry, ALPHA specifies the scalar alpha. +* Unchanged on exit. +* +* A - DOUBLE PRECISION array of DIMENSION ( LDA, ka ), where ka is +* k when TRANSA = 'N' or 'n', and is m otherwise. +* Before entry with TRANSA = 'N' or 'n', the leading m by k +* part of the array A must contain the matrix A, otherwise +* the leading k by m part of the array A must contain the +* matrix A. +* Unchanged on exit. +* +* LDA - INTEGER. +* On entry, LDA specifies the first dimension of A as declared +* in the calling (sub) program. When TRANSA = 'N' or 'n' then +* LDA must be at least max( 1, m ), otherwise LDA must be at +* least max( 1, k ). +* Unchanged on exit. +* +* B - DOUBLE PRECISION array of DIMENSION ( LDB, kb ), where kb is +* n when TRANSB = 'N' or 'n', and is k otherwise. +* Before entry with TRANSB = 'N' or 'n', the leading k by n +* part of the array B must contain the matrix B, otherwise +* the leading n by k part of the array B must contain the +* matrix B. +* Unchanged on exit. +* +* LDB - INTEGER. +* On entry, LDB specifies the first dimension of B as declared +* in the calling (sub) program. When TRANSB = 'N' or 'n' then +* LDB must be at least max( 1, k ), otherwise LDB must be at +* least max( 1, n ). +* Unchanged on exit. +* +* BETA - DOUBLE PRECISION. +* On entry, BETA specifies the scalar beta. When BETA is +* supplied as zero then C need not be set on input. +* Unchanged on exit. +* +* C - DOUBLE PRECISION array of DIMENSION ( LDC, n ). +* Before entry, the leading m by n part of the array C must +* contain the matrix C, except when beta is zero, in which +* case C need not be set on entry. +* On exit, the array C is overwritten by the m by n matrix +* ( alpha*op( A )*op( B ) + beta*C ). +* +* LDC - INTEGER. +* On entry, LDC specifies the first dimension of C as declared +* in the calling (sub) program. LDC must be at least +* max( 1, m ). +* Unchanged on exit. +* +* +* Level 3 Blas routine. +* +* -- Written on 8-February-1989. +* Jack Dongarra, Argonne National Laboratory. +* Iain Duff, AERE Harwell. +* Jeremy Du Croz, Numerical Algorithms Group Ltd. +* Sven Hammarling, Numerical Algorithms Group Ltd. +* +* +* .. External Functions .. + LOGICAL LSAME + EXTERNAL LSAME +* .. External Subroutines .. + EXTERNAL XERBLA +* .. Intrinsic Functions .. + INTRINSIC MAX +* .. Local Scalars .. + LOGICAL NOTA, NOTB + INTEGER I, INFO, J, L, NCOLA, NROWA, NROWB + DOUBLE PRECISION TEMP +* .. Parameters .. + DOUBLE PRECISION ONE , ZERO + PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) +* .. +* .. Executable Statements .. +* +* Set NOTA and NOTB as true if A and B respectively are not +* transposed and set NROWA, NCOLA and NROWB as the number of rows +* and columns of A and the number of rows of B respectively. +* + NOTA = LSAME( TRANSA, 'N' ) + NOTB = LSAME( TRANSB, 'N' ) + IF( NOTA )THEN + NROWA = M + NCOLA = K + ELSE + NROWA = K + NCOLA = M + END IF + IF( NOTB )THEN + NROWB = K + ELSE + NROWB = N + END IF +* +* Test the input parameters. +* + INFO = 0 + IF( ( .NOT.NOTA ).AND. + $ ( .NOT.LSAME( TRANSA, 'C' ) ).AND. + $ ( .NOT.LSAME( TRANSA, 'T' ) ) )THEN + INFO = 1 + ELSE IF( ( .NOT.NOTB ).AND. + $ ( .NOT.LSAME( TRANSB, 'C' ) ).AND. + $ ( .NOT.LSAME( TRANSB, 'T' ) ) )THEN + INFO = 2 + ELSE IF( M .LT.0 )THEN + INFO = 3 + ELSE IF( N .LT.0 )THEN + INFO = 4 + ELSE IF( K .LT.0 )THEN + INFO = 5 + ELSE IF( LDA.LT.MAX( 1, NROWA ) )THEN + INFO = 8 + ELSE IF( LDB.LT.MAX( 1, NROWB ) )THEN + INFO = 10 + ELSE IF( LDC.LT.MAX( 1, M ) )THEN + INFO = 13 + END IF + IF( INFO.NE.0 )THEN + CALL XERBLA( 'DGEMM ', INFO ) + RETURN + END IF +* +* Quick return if possible. +* + IF( ( M.EQ.0 ).OR.( N.EQ.0 ).OR. + $ ( ( ( ALPHA.EQ.ZERO ).OR.( K.EQ.0 ) ).AND.( BETA.EQ.ONE ) ) ) + $ RETURN +* +* And if alpha.eq.zero. +* + IF( ALPHA.EQ.ZERO )THEN + IF( BETA.EQ.ZERO )THEN + DO 20, J = 1, N + DO 10, I = 1, M + C( I, J ) = ZERO + 10 CONTINUE + 20 CONTINUE + ELSE + DO 40, J = 1, N + DO 30, I = 1, M + C( I, J ) = BETA*C( I, J ) + 30 CONTINUE + 40 CONTINUE + END IF + RETURN + END IF +* +* Start the operations. +* + IF( NOTB )THEN + IF( NOTA )THEN +* +* Form C := alpha*A*B + beta*C. +* + DO 90, J = 1, N + IF( BETA.EQ.ZERO )THEN + DO 50, I = 1, M + C( I, J ) = ZERO + 50 CONTINUE + ELSE IF( BETA.NE.ONE )THEN + DO 60, I = 1, M + C( I, J ) = BETA*C( I, J ) + 60 CONTINUE + END IF + DO 80, L = 1, K + IF( B( L, J ).NE.ZERO )THEN + TEMP = ALPHA*B( L, J ) + DO 70, I = 1, M + C( I, J ) = C( I, J ) + TEMP*A( I, L ) + 70 CONTINUE + END IF + 80 CONTINUE + 90 CONTINUE + ELSE +* +* Form C := alpha*A'*B + beta*C +* + DO 120, J = 1, N + DO 110, I = 1, M + TEMP = ZERO + DO 100, L = 1, K + TEMP = TEMP + A( L, I )*B( L, J ) + 100 CONTINUE + IF( BETA.EQ.ZERO )THEN + C( I, J ) = ALPHA*TEMP + ELSE + C( I, J ) = ALPHA*TEMP + BETA*C( I, J ) + END IF + 110 CONTINUE + 120 CONTINUE + END IF + ELSE + IF( NOTA )THEN +* +* Form C := alpha*A*B' + beta*C +* + DO 170, J = 1, N + IF( BETA.EQ.ZERO )THEN + DO 130, I = 1, M + C( I, J ) = ZERO + 130 CONTINUE + ELSE IF( BETA.NE.ONE )THEN + DO 140, I = 1, M + C( I, J ) = BETA*C( I, J ) + 140 CONTINUE + END IF + DO 160, L = 1, K + IF( B( J, L ).NE.ZERO )THEN + TEMP = ALPHA*B( J, L ) + DO 150, I = 1, M + C( I, J ) = C( I, J ) + TEMP*A( I, L ) + 150 CONTINUE + END IF + 160 CONTINUE + 170 CONTINUE + ELSE +* +* Form C := alpha*A'*B' + beta*C +* + DO 200, J = 1, N + DO 190, I = 1, M + TEMP = ZERO + DO 180, L = 1, K + TEMP = TEMP + A( L, I )*B( J, L ) + 180 CONTINUE + IF( BETA.EQ.ZERO )THEN + C( I, J ) = ALPHA*TEMP + ELSE + C( I, J ) = ALPHA*TEMP + BETA*C( I, J ) + END IF + 190 CONTINUE + 200 CONTINUE + END IF + END IF +* + RETURN +* +* End of DGEMM . +* + END + SUBROUTINE DGEMV ( TRANS, M, N, ALPHA, A, LDA, X, INCX, + $ BETA, Y, INCY ) +* .. Scalar Arguments .. + DOUBLE PRECISION ALPHA, BETA + INTEGER INCX, INCY, LDA, M, N + CHARACTER*1 TRANS +* .. Array Arguments .. + DOUBLE PRECISION A( LDA, * ), X( * ), Y( * ) +* .. +* +* Purpose +* ======= +* +* DGEMV performs one of the matrix-vector operations +* +* y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y, +* +* where alpha and beta are scalars, x and y are vectors and A is an +* m by n matrix. +* +* Parameters +* ========== +* +* TRANS - CHARACTER*1. +* On entry, TRANS specifies the operation to be performed as +* follows: +* +* TRANS = 'N' or 'n' y := alpha*A*x + beta*y. +* +* TRANS = 'T' or 't' y := alpha*A'*x + beta*y. +* +* TRANS = 'C' or 'c' y := alpha*A'*x + beta*y. +* +* Unchanged on exit. +* +* M - INTEGER. +* On entry, M specifies the number of rows of the matrix A. +* M must be at least zero. +* Unchanged on exit. +* +* N - INTEGER. +* On entry, N specifies the number of columns of the matrix A. +* N must be at least zero. +* Unchanged on exit. +* +* ALPHA - DOUBLE PRECISION. +* On entry, ALPHA specifies the scalar alpha. +* Unchanged on exit. +* +* A - DOUBLE PRECISION array of DIMENSION ( LDA, n ). +* Before entry, the leading m by n part of the array A must +* contain the matrix of coefficients. +* Unchanged on exit. +* +* LDA - INTEGER. +* On entry, LDA specifies the first dimension of A as declared +* in the calling (sub) program. LDA must be at least +* max( 1, m ). +* Unchanged on exit. +* +* X - DOUBLE PRECISION array of DIMENSION at least +* ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n' +* and at least +* ( 1 + ( m - 1 )*abs( INCX ) ) otherwise. +* Before entry, the incremented array X must contain the +* vector x. +* Unchanged on exit. +* +* INCX - INTEGER. +* On entry, INCX specifies the increment for the elements of +* X. INCX must not be zero. +* Unchanged on exit. +* +* BETA - DOUBLE PRECISION. +* On entry, BETA specifies the scalar beta. When BETA is +* supplied as zero then Y need not be set on input. +* Unchanged on exit. +* +* Y - DOUBLE PRECISION array of DIMENSION at least +* ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n' +* and at least +* ( 1 + ( n - 1 )*abs( INCY ) ) otherwise. +* Before entry with BETA non-zero, the incremented array Y +* must contain the vector y. On exit, Y is overwritten by the +* updated vector y. +* +* INCY - INTEGER. +* On entry, INCY specifies the increment for the elements of +* Y. INCY must not be zero. +* Unchanged on exit. +* +* +* Level 2 Blas routine. +* +* -- Written on 22-October-1986. +* Jack Dongarra, Argonne National Lab. +* Jeremy Du Croz, Nag Central Office. +* Sven Hammarling, Nag Central Office. +* Richard Hanson, Sandia National Labs. +* +* +* .. Parameters .. + DOUBLE PRECISION ONE , ZERO + PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) +* .. Local Scalars .. + DOUBLE PRECISION TEMP + INTEGER I, INFO, IX, IY, J, JX, JY, KX, KY, LENX, LENY +* .. External Functions .. + LOGICAL LSAME + EXTERNAL LSAME +* .. External Subroutines .. + EXTERNAL XERBLA +* .. Intrinsic Functions .. + INTRINSIC MAX +* .. +* .. Executable Statements .. +* +* Test the input parameters. +* + INFO = 0 + IF ( .NOT.LSAME( TRANS, 'N' ).AND. + $ .NOT.LSAME( TRANS, 'T' ).AND. + $ .NOT.LSAME( TRANS, 'C' ) )THEN + INFO = 1 + ELSE IF( M.LT.0 )THEN + INFO = 2 + ELSE IF( N.LT.0 )THEN + INFO = 3 + ELSE IF( LDA.LT.MAX( 1, M ) )THEN + INFO = 6 + ELSE IF( INCX.EQ.0 )THEN + INFO = 8 + ELSE IF( INCY.EQ.0 )THEN + INFO = 11 + END IF + IF( INFO.NE.0 )THEN + CALL XERBLA( 'DGEMV ', INFO ) + RETURN + END IF +* +* Quick return if possible. +* + IF( ( M.EQ.0 ).OR.( N.EQ.0 ).OR. + $ ( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) ) + $ RETURN +* +* Set LENX and LENY, the lengths of the vectors x and y, and set +* up the start points in X and Y. +* + IF( LSAME( TRANS, 'N' ) )THEN + LENX = N + LENY = M + ELSE + LENX = M + LENY = N + END IF + IF( INCX.GT.0 )THEN + KX = 1 + ELSE + KX = 1 - ( LENX - 1 )*INCX + END IF + IF( INCY.GT.0 )THEN + KY = 1 + ELSE + KY = 1 - ( LENY - 1 )*INCY + END IF +* +* Start the operations. In this version the elements of A are +* accessed sequentially with one pass through A. +* +* First form y := beta*y. +* + IF( BETA.NE.ONE )THEN + IF( INCY.EQ.1 )THEN + IF( BETA.EQ.ZERO )THEN + DO 10, I = 1, LENY + Y( I ) = ZERO + 10 CONTINUE + ELSE + DO 20, I = 1, LENY + Y( I ) = BETA*Y( I ) + 20 CONTINUE + END IF + ELSE + IY = KY + IF( BETA.EQ.ZERO )THEN + DO 30, I = 1, LENY + Y( IY ) = ZERO + IY = IY + INCY + 30 CONTINUE + ELSE + DO 40, I = 1, LENY + Y( IY ) = BETA*Y( IY ) + IY = IY + INCY + 40 CONTINUE + END IF + END IF + END IF + IF( ALPHA.EQ.ZERO ) + $ RETURN + IF( LSAME( TRANS, 'N' ) )THEN +* +* Form y := alpha*A*x + y. +* + JX = KX + IF( INCY.EQ.1 )THEN + DO 60, J = 1, N + IF( X( JX ).NE.ZERO )THEN + TEMP = ALPHA*X( JX ) + DO 50, I = 1, M + Y( I ) = Y( I ) + TEMP*A( I, J ) + 50 CONTINUE + END IF + JX = JX + INCX + 60 CONTINUE + ELSE + DO 80, J = 1, N + IF( X( JX ).NE.ZERO )THEN + TEMP = ALPHA*X( JX ) + IY = KY + DO 70, I = 1, M + Y( IY ) = Y( IY ) + TEMP*A( I, J ) + IY = IY + INCY + 70 CONTINUE + END IF + JX = JX + INCX + 80 CONTINUE + END IF + ELSE +* +* Form y := alpha*A'*x + y. +* + JY = KY + IF( INCX.EQ.1 )THEN + DO 100, J = 1, N + TEMP = ZERO + DO 90, I = 1, M + TEMP = TEMP + A( I, J )*X( I ) + 90 CONTINUE + Y( JY ) = Y( JY ) + ALPHA*TEMP + JY = JY + INCY + 100 CONTINUE + ELSE + DO 120, J = 1, N + TEMP = ZERO + IX = KX + DO 110, I = 1, M + TEMP = TEMP + A( I, J )*X( IX ) + IX = IX + INCX + 110 CONTINUE + Y( JY ) = Y( JY ) + ALPHA*TEMP + JY = JY + INCY + 120 CONTINUE + END IF + END IF +* + RETURN +* +* End of DGEMV . +* + END + + subroutine perma(compo) + include "calcium.hf" + double precision a(5,24),b(24),q(6) + dimension t(24),tn(24) + dimension puissn(24),puiss(24) + dimension text(6),rflu(6),textn(6) + dimension tpar(6),rpar(6),tparn(6) + integer compo + + ldab=5 + ldb=24 + nrhs=1 + npt=0 + r=1. + ro=1. + rext=0.5 + te=10. + do i=1,24 + tn(i)=100. + t(i)=100. + puiss(i)=100. + enddo + do j=1,6 + q(j)=50. + enddo +c +c construction de la matrice (laplacien) +c + do i = 1, 5 + do j = 1, 24 + a(i,j) = 0. + enddo + enddo + + do i=2,3 + do j=2,5 + npt=i+(j-1)*4 + a(1,npt)=4./r + a(2,npt)=-1./r + a(5,npt)=-1./r + enddo + enddo + do i=2,3 + npt=i + a(1,npt)=3./r + a(2,npt)=-1./r + a(5,npt)=-1./r + npt=i+20 + a(1,npt)=3./r + a(2,npt)=-1./r + a(5,npt)=0. + enddo + do j=2,5 + npt=1+(j-1)*4 + a(1,npt)=3./r + a(2,npt)=-1./r + a(5,npt)=-1./r + npt=4+(j-1)*4 + a(1,npt)=3./r+1./(r/2.+rext) + a(2,npt)=0. + a(5,npt)=-1./r + enddo + i=1 + a(1,i)=2./r + a(2,i)=-1./r + a(5,i)=-1./r + i=21 + a(1,i)=2./r + a(2,i)=-1./r + a(5,i)=-1./r + i=24 + a(1,i)=2./r+1./(r/2.+rext) + a(2,i)=-1./r + a(5,i)=-1./r + i=4 + a(1,i)=2./r+1./(r/2.+rext) + a(2,i)=0. + a(5,i)=-1./r + +c +c factorisation de la matrice +c + n=24 + kd=4 + + call dPBTRF( 'L' , N, KD, A, LDAB, INFO ) + + iter=0 + + do i=1,6 + tpar(i)=t(4*i) + rpar(i)=r/2. + enddo + + do i=1,24 + tn(i)=0. + puissn(i)=0. + enddo + do i=1,6 + tparn(i)=0. + textn(i)=0. + rflu(i)=0. + enddo +c +c initialisation de la temperature a iter=0 +c + CALL cpeRE(compo,CP_ITERATION, ti, iter, 'tempi', 24, + & t , info) + IF( info.EQ.CPSTOP )GO TO 9000 +c +c initialisation de la temperature de bord a iter=0. +c + CALL cpeRE(compo,CP_ITERATION, ti, iter, 'tpi', 6, + & tpar , info) + IF( info.EQ.CPSTOP )GO TO 9000 +c +c boucle iterative infinie +c + do while( iter .lt. 500) +c +c lecture de la puissance combustible a iter +c + CALL cplRE(compo,CP_ITERATION,ti, tf, iter, 'puissi', 24, + & nval, puiss , info) + IF( info.EQ.CPSTOP )GO TO 9000 +c +c lecture de la temperature exterieure a iter +c + CALL cplRE(compo,CP_ITERATION,ti, tf, iter, 'tfi', 6, + & nval, text , info) + IF( info.EQ.CPSTOP )GO TO 9000 +c +c calcul du second membre +c + do npt=1,24 + b(npt)=puiss(npt) + enddo + + do j=1,6 + npt=j*4 + b(npt)=b(npt)+text(j)/(r/2+rflu(j)) + enddo + +c +c resolution du systeme lineaire +c + call dPBTRs( 'L' , N,kd, nrhs,A,LDAB,b,ldb,INFO ) + + do i=1,24 + t(i)=b(i) + enddo + + do i=1,6 + tpar(i)=t(4*i) + enddo + + iter=iter+1 +c +c ecriture de la temperature a iter+1 +c + CALL cpeRE(compo,CP_ITERATION, ti, iter, 'tempi', 24, + & t , info) + IF( info.NE. CPOK )GO TO 9000 +c +c ecriture de la temperature de paroi a iter+1 +c + CALL cpeRE(compo,CP_ITERATION, ti, iter, 'tpi', 6, + & tpar , info) + IF( info.NE. CPOK )GO TO 9000 +c +c calcul de convergence +c + conv1=0. + do i=1,npt + conv1=conv1+(puiss(i)-puissn(i))**2 + conv1=conv1+(t(i)-tn(i))**2 + enddo + do i=1,6 + conv1=conv1+(tpar(i)-tparn(i))**2 + conv1=conv1+(text(i)-textn(i))**2 + enddo + + iconv=0 + if(conv1.lt.1.e-3) iconv=1 +c +c ecriture du flag de convergence iconv +c + write(6,*)"SOLIDE:",iter,iconv + call flush(6) + CALL cpeEN(compo,CP_ITERATION,tf, iter , 'iconv', 1, + & iconv , info) + IF( info.NE. CPOK )GO TO 9000 + + if(iconv.eq.1)go to 9000 + + do i=1,24 + tn(i)=b(i) + puissn(i)=puiss(i) + enddo + do i=1,6 + tparn(i)=tpar(i) + textn(i)=text(i) + enddo + + enddo + +9000 continue + CALL cpfin(compo,CP_ARRET, info) + end -- 2.30.2