scnc
/
qmcpack
mirror of https://github.com/QMCPACK/qmcpack.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Merge branch 'develop' into UnifiedDriver_branch_engine

This commit is contained in:
Peter Doak 2019-10-01 14:28:01 -04:00
parent 40a586503f b1ba44ece4
commit 007bf641cb
47 changed files with 2473 additions and 1859 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

5
external_codes/mpi_wrapper/mpi3/communicator.hpp
View File

@ -17,6 +17,7 @@
#include "../mpi3/status.hpp"
#include "../mpi3/type.hpp"
#include "../mpi3/error.hpp"
#include "../mpi3/group.hpp"
#include "../mpi3/detail/basic_communicator.hpp"
#include "../mpi3/detail/buffer.hpp"
@ -146,8 +147,6 @@ class group;
using boost::optional;
template<class T = void>
struct window;
struct error_handler;
struct keyval;
@ -185,6 +184,8 @@ struct shared_communicator; // intracommunicator
using boost::any;
using boost::any_cast;
template<class T> struct window;
class communicator : protected detail::basic_communicator{
friend struct detail::package;
protected:

2
external_codes/mpi_wrapper/mpi3/dynamic_window.hpp
View File

@ -16,7 +16,7 @@ struct dynamic_window : window<T>{
dynamic_window() : window<T>{}{}
public:
dynamic_window(communicator& comm){
int s = MPI_Win_create_dynamic(MPI_INFO_NULL, comm.impl_, &(this->impl_));
int s = MPI_Win_create_dynamic(MPI_INFO_NULL, comm.get(), &(this->impl_));
if(s != MPI_SUCCESS) throw std::runtime_error("cannot create dynamic window");
}
template<class TT = char>

69
external_codes/mpi_wrapper/mpi3/shm/allocator.hpp
View File

@ -4,22 +4,57 @@
#ifndef MPI3_SHM_ALLOCATOR_HPP
#define MPI3_SHM_ALLOCATOR_HPP
#include "../../mpi3/shared_window.hpp"
//#include "../../mpi3/shared_window.hpp"
#include "../shm/pointer.hpp"
namespace boost{
namespace mpi3{
namespace shm{
template<class... As>
using allocator = mpi3::intranode::allocator<As...>;
//template<class... As>
//using allocator = mpi3::intranode::allocator<As...>;
template<class T>
using pointer = mpi3::intranode::array_ptr<T>;
//template<class T>
//using pointer = mpi3::intranode::array_ptr<T>;
template<class Ptr>
struct pointer_traits : std::pointer_traits<Ptr>{
static auto to_address(Ptr const& p){
return std::addressof(*p);
template<class T = void>
struct allocator{
template<class U> struct rebind{typedef allocator<U> other;};
using value_type = T;
using pointer = shm::pointer<T>;
using const_pointer = shm::pointer<T const>;
using size_type = mpi3::size_t;
using difference_type = mpi3::ptrdiff_t;
private:
mpi3::shared_communicator* comm_;
allocator() = delete;
public:
allocator(mpi3::shared_communicator* comm) : comm_{comm}{}
template<class U> allocator(allocator<U> const& o) : comm_{o.comm_}{}
~allocator() = default;
pointer allocate(size_type n, const void* /*hint*/ = 0) const{
pointer ret = 0;
// this code using reset is to workaround a bug in outdated libc++ https://bugs.llvm.org/show_bug.cgi?id=42021
ret.w_.reset(new shared_window<T>{
comm_->make_shared_window<T>(comm_->root()?n:0)
});
// ret.w_ = std::make_shared<shared_window<T>>(
// comm_->make_shared_window<T>(comm_->root()?n:0)
// );
return ret;
}
void deallocate(pointer p, size_type){p.w_.reset();}
allocator& operator=(allocator const& other) = default;
bool operator==(allocator const& other) const{return comm_==other.comm_;}
bool operator!=(allocator const& other) const{return not(other == *this);}
template<class P = pointer, class... As>
void construct(P p, As&&... as){
if(comm_->root()) ::new((void*)p) typename std::pointer_traits<P>::element_type(std::forward<As>(as)...);
}
template<class P = pointer> void destroy(P p){
using V = typename std::pointer_traits<P>::element_type;
p->~V();
}
};
@ -33,6 +68,7 @@ struct pointer_traits : std::pointer_traits<Ptr>{
namespace mpi3 = boost::mpi3;
#if 0
template<class A>
struct monotonic{
A& a_;
@ -57,22 +93,13 @@ struct monotonic{
}
void deallocate(monotonic::pointer, monotonic::size_type){}
};
#endif
int mpi3::main(int, char*[], mpi3::communicator world){
std::allocator<double> stda;
std::allocator<double>::rebind<int>::other dd;
monotonic<std::allocator<double>> a(1024, stda);
std::vector<double, std::allocator<double>> v1(600, stda);
std::vector<double, monotonic<std::allocator<double> > > v2(600, a);
#if 1
// std::vector<double, monotonic<std::allocator<double>>> v2(600, a);
mpi3::shared_communicator node = world.split_shared();
mpi3::shm::allocator<double> A1(node);
mpi3::shm::allocator<double> A1(std::addressof(node));
mpi3::shm::pointer<double> data1 = A1.allocate(80);
using ptr = decltype(data1);
@ -90,7 +117,7 @@ int mpi3::main(int, char*[], mpi3::communicator world){
assert(data1[3] == 3.4);
A1.deallocate(data1, 80);
#endif
return 0;
}

17
external_codes/mpi_wrapper/mpi3/shm/memory.hpp
View File

@ -81,6 +81,23 @@ struct mapped_region{
}}}
namespace boost{
namespace mpi3{
namespace shm{
template<class Alloc, typename T, typename Size, typename TT>
pointer<T> uninitialized_fill_n(Alloc const& a, pointer<T> f, Size n, TT const& val){
assert(0);
using std::uninitialized_fill_n;
if(mpi3::group(*f.wSP_).root()) uninitialized_fill_n(to_address(f), n, val);
first.wSP_->fence();
first.wSP_->fence();
return first + n;
}
}}}
#ifdef _TEST_BOOST_MPI3_SHM_MEMORY
#include "../../mpi3/main.hpp"

6
external_codes/mpi_wrapper/mpi3/shm/mutex.hpp
View File

@ -15,14 +15,14 @@ class mutex{
mpi3::shm::allocator<std::atomic_flag> alloc_;//(node);
mpi3::shm::pointer<std::atomic_flag> f_;
public:
mutex(mpi3::shared_communicator& scomm) : scomm_(scomm), alloc_(scomm_), f_(alloc_.allocate(1)){
mutex(mpi3::shared_communicator& scomm) : scomm_(scomm), alloc_(std::addressof(scomm_)), f_(alloc_.allocate(1)){
if(scomm_.root()) alloc_.construct(&*f_, false);
scomm_.barrier();
}
mutex(mutex const&) = delete;
mutex(mutex&&) = delete;
void lock(){while(f_->test_and_set(std::memory_order_acquire));}
void unlock(){f_->clear(std::memory_order_release);}
void lock(){while((*f_).test_and_set(std::memory_order_acquire));}
void unlock(){(*f_).clear(std::memory_order_release);}
~mutex(){
if(scomm_.root()) alloc_.destroy(&*f_);
scomm_.barrier();

181
external_codes/mpi_wrapper/mpi3/window.hpp
View File

@ -1,12 +1,15 @@
#if COMPILATION_INSTRUCTIONS
(echo "#include\""$0"\"" > $0x.cpp) && mpic++ -O3 -std=c++14 -Wfatal-errors -D_TEST_BOOST_MPI3_WINDOW $0x.cpp -o $0x.x && time mpirun -np 4 $0x.x $@ && rm -f $0x.x $0x.cpp; exit
(echo "#include\""$0"\"" > $0x.cpp) && mpic++ -O3 -std=c++17 `#-Wfatal-errors` -Wall -Wextra -Wpedantic -D_TEST_BOOST_MPI3_WINDOW $0x.cpp -o $0x.x && time mpirun -np 4 $0x.x $@ && rm -f $0x.x $0x.cpp; exit
#endif
// (C) Copyright Alfredo A. Correa 2018.
// (C) Copyright Alfredo A. Correa 2019
#ifndef BOOST_MPI3_WINDOW_HPP
#define BOOST_MPI3_WINDOW_HPP
#include "../mpi3/communicator.hpp"
#include "../mpi3/error.hpp"
#include "../mpi3/type.hpp"
#include "../mpi3/detail/datatype.hpp"
#include "../mpi3/detail/call.hpp"
#include "../mpi3/communicator.hpp"
#define OMPI_SKIP_MPICXX 1 // https://github.com/open-mpi/ompi/issues/5157
#include<mpi.h>
@ -14,51 +17,42 @@
namespace boost{
namespace mpi3{
class group;
struct target{
int rank;
mpi3::size_t disp;
};
template<class T = void>
class panel;
template<class T = void> class panel;
//template<class T = void> struct window;
template<class T = void> struct window;
struct basic_window{
template<>
struct window<void>{
protected:
MPI_Win impl_ = MPI_WIN_NULL;
basic_window() = default;
basic_window(MPI_Win w) : impl_{w}{}
public:
MPI_Win& operator&(){return impl_;}
MPI_Win const& operator&() const{return impl_;}
};
template<>
struct window<void> : basic_window{
public:
window(){}
template<class T>
window(T* base, mpi3::size_t size, communicator& comm){
int s = MPI_Win_create(
(void*)base, size*sizeof(T), sizeof(T), MPI_INFO_NULL,
&comm, &impl_
); // TODO probably need alignof
if(s != MPI_SUCCESS) throw std::runtime_error{"cannot create window"};
void clear(){
if(impl_ != MPI_WIN_NULL) MPI3_CALL(MPI_Win_free)(&impl_);
assert(impl_ == MPI_WIN_NULL);
}
window(void* base, mpi3::size_t size, communicator& comm){
int s = MPI_Win_create(base, size, 1, MPI_INFO_NULL, &comm, &impl_);
if(s != MPI_SUCCESS) throw std::runtime_error{"cannot create window"};
protected:
window() = default;
public:
template<class T, class Size = mpi3::size_t>
window(communicator const& c, T* b, Size n = 0){
MPI3_CALL(MPI_Win_create)(b, n*sizeof(T), alignof(T), MPI_INFO_NULL, c.get(), &impl_);
assert( alignof(T) == sizeof(T) ); // to see in what type it is different
}
window(communicator& comm) : window((void*)nullptr, 0, comm){}
window(window const&) = delete; // cannot be duplicated, see text before sec. 4.5 in Using Adv. MPI
window(window&& other) : basic_window{std::exchange(other.impl_, MPI_WIN_NULL)}{//is movable if null is not a correct state?
// other.impl_ = MPI_WIN_NULL;
window(window const&) = delete;// see text before §4.5 in Using Adv. MPI
window(window&& o) noexcept : impl_{std::exchange(o.impl_, MPI_WIN_NULL)}{}
window& operator=(window const&) = delete; // see cctor
window& operator=(window&& other){// self assignment is undefined
clear(); swap(*this, other); return *this;
}
window& operator=(window const&) = delete; // see comment in cctor
window& operator=(window&& other){ // guard self assignment? probably worthless
if(impl_ != MPI_WIN_NULL) MPI_Win_free(&impl_);
impl_ = std::exchange(other.impl_, MPI_WIN_NULL);
return *this;
}
~window(){if(impl_ != MPI_WIN_NULL) MPI_Win_free(&impl_);}
friend void swap(window& a, window& b){std::swap(a.impl_, b.impl_);}
~window(){clear();}
template<typename It1, typename Size, typename V = typename std::iterator_traits<It1>::value_type>
void accumulate_n(It1 first, Size count, int target_rank, int target_disp = 0){
using detail::data;
@ -84,37 +78,40 @@ struct window<void> : basic_window{
void flush_local(){return flush_local_all();}
void* base() const{
void* base; int flag;
int s = MPI_Win_get_attr(impl_, MPI_WIN_BASE, &base, &flag);
if(s != MPI_SUCCESS) throw std::runtime_error("cannot get base");
MPI3_CALL(MPI_Win_get_attr)(impl_, MPI_WIN_BASE, &base, &flag);
assert(flag);
return base;
}
mpi3::size_t const& size() const{
MPI_Aint* size_p; int flag;
int s = MPI_Win_get_attr(impl_, MPI_WIN_SIZE, &size_p, &flag);
if(s != MPI_SUCCESS) throw std::runtime_error("cannot get base");
MPI3_CALL(MPI_Win_get_attr)(impl_, MPI_WIN_SIZE, &size_p, &flag);
assert(flag);
return *size_p;
}
int const& disp_unit() const{
int* disp_unit_p; int flag;
int s = MPI_Win_get_attr(impl_, MPI_WIN_DISP_UNIT, &disp_unit_p, &flag);
if(s != MPI_SUCCESS) throw std::runtime_error("cannot get base");
MPI3_CALL(MPI_Win_get_attr)(impl_, MPI_WIN_DISP_UNIT, &disp_unit_p, &flag);
assert(flag);
return *disp_unit_p;
}
// get_errhandler(...);
group get_group() const{
group ret;
MPI3_CALL(MPI_Win_get_group)(impl_, &(&ret));
return ret;
}
// group get_group(){use reinterpret_cast?}
// ... get_info
// ... get_name
// lock arguments are reversed
void lock(int rank, int lock_type = MPI_LOCK_EXCLUSIVE, int assert = MPI_MODE_NOCHECK){
MPI_Win_lock(lock_type, rank, assert, impl_);
MPI3_CALL(MPI_Win_lock)(lock_type, rank, assert, impl_);
}
void lock_exclusive(int rank, int assert = MPI_MODE_NOCHECK){
MPI_Win_lock(MPI_LOCK_EXCLUSIVE, rank, assert, impl_);
MPI3_CALL(MPI_Win_lock)(MPI_LOCK_EXCLUSIVE, rank, assert, impl_);
}
void lock_shared(int rank, int assert = MPI_MODE_NOCHECK){
MPI_Win_lock(MPI_LOCK_SHARED, rank, assert, impl_);
MPI3_CALL(MPI_Win_lock)(MPI_LOCK_SHARED, rank, assert, impl_);
}
void lock_all(int assert = MPI_MODE_NOCHECK){MPI_Win_lock_all(assert, impl_);}
// void post(group const& g, int assert = MPI_MODE_NOCHECK) const{
@ -133,39 +130,33 @@ struct window<void> : basic_window{
void unlock(int rank) const{MPI_Win_unlock(rank, impl_);}
void unlock_all(){MPI_Win_unlock_all(impl_);}
void wait() const{MPI_Win_wait(impl_);}
// void fetch_and_op(T const* origin, T* target, int target_rank, int target_disp = 0) const{
// MPI_Fetch_and_op(origin, target, detail::datatype<T>{}, target_rank, target_disp, , impl_);
// }
// template<class T, class Op, class datatype = detail::datatype<T>, >
// void fetch_and_op(T const* origin, T* target, int target_rank, int target_disp = 0) const{
// MPI_Fetch_and_op(origin, target, datatype{}, target_rank, target_disp, , impl_);
// }
// void fetch_exchange(T const* origin, T* target, int target_rank, int target_disp = 0) const{
// MPI_Fetch_and_op(origin, target,detail::datatype<T>{}, target_rank, target_disp, MPI_REPLACE, impl_);
// }
// maybe this goes to a pointer impl
template<class T>
void fetch_sum_value(T const& origin, T& target, int target_rank, int target_disp = 0) const{
MPI_Fetch_and_op(&origin, &target, detail::basic_datatype<T>{}, target_rank, target_disp, MPI_SUM, impl_);
void fetch_sum_value(T const& origin, T& target, int target_rank, int target_disp=0) const{
MPI3_CALL(MPI_Fetch_and_op)(&origin, &target, detail::basic_datatype<T>{}, target_rank, target_disp, MPI_SUM, impl_);
}
template<class T>
void fetch_prod_value(T const& origin, T& target, int target_rank, int target_disp = 0) const{
MPI_Fetch_and_op(&origin, &target, detail::basic_datatype<T>{}, target_rank, target_disp, MPI_PROD, impl_);
MPI3_CALL(MPI_Fetch_and_op)(&origin, &target, detail::basic_datatype<T>{}, target_rank, target_disp, MPI_PROD, impl_);
}
template<class T>
void fetch_replace_value(T const& origin, T& target, int target_rank, int target_disp = 0) const{
MPI_Fetch_and_op(&origin, &target, detail::basic_datatype<T>{}, target_rank, target_disp, MPI_REPLACE, impl_);
MPI3_CALL(MPI_Fetch_and_op)(&origin, &target, detail::basic_datatype<T>{}, target_rank, target_disp, MPI_REPLACE, impl_);
}
template<class CI1, class CI2, class datatype = detail::basic_datatype<typename std::iterator_traits<CI1>::value_type> >
void fetch_replace(CI1 it1, CI2 it2, int target_rank, int target_disp = 0) const{
MPI_Fetch_and_op(std::addressof(*it1), std::addressof(*it2), datatype{}, target_rank, target_disp, MPI_REPLACE, impl_);
MPI3_CALL(MPI_Fetch_and_op)(std::addressof(*it1), std::addressof(*it2), datatype{}, target_rank, target_disp, MPI_REPLACE, impl_);
}
template<class ContiguousIterator>
void blocking_put_n(ContiguousIterator it, int count, int target_rank, int target_offset = 0){
lock_shared(target_rank, 0);
@ -175,7 +166,7 @@ struct window<void> : basic_window{
template<class ContiguousIterator>
void put_n(ContiguousIterator it, std::size_t n, int target_rank, int target_disp = 0) const{
using detail::data;
MPI_Put(
MPI3_CALL(MPI_Put)(
data(it), /* void* origin_address = a + i*/
n, /*int origin_count = 1 */
detail::basic_datatype<typename std::iterator_traits<ContiguousIterator>::value_type>{},
@ -198,7 +189,7 @@ struct window<void> : basic_window{
template<typename ContiguousIterator, typename Size>
ContiguousIterator get_n(ContiguousIterator it, Size n, int target_rank, int target_disp = 0) const{
using detail::data;
int s = MPI_Get(
MPI3_CALL(MPI_Get)(
data(it), /* void* origin_address = b + i*/
n, /*int origin_count = 1 */
detail::basic_datatype<typename std::iterator_traits<ContiguousIterator>::value_type>{},
@ -208,7 +199,6 @@ struct window<void> : basic_window{
detail::basic_datatype<typename std::iterator_traits<ContiguousIterator>::value_type>{},
impl_
);
if(s != MPI_SUCCESS) throw std::runtime_error("cannot get_n");
return it + n;
}
template<typename ContiguousIterator>
@ -224,9 +214,12 @@ struct window<void> : basic_window{
template<class T>
struct window : window<void>{
protected:
window() = default;
window(T* base, mpi3::size_t size, communicator& comm) : window<void>(static_cast<void*>(base), size*sizeof(T), comm){}
T* base() const{return window<void>::base();}
public:
template<class Size = mpi3::size_t>
window(communicator const& c, T* b, Size n = 0) : window<void>{c, b, n}{}
T* base() const{return static_cast<T*>(window<void>::base());}
mpi3::size_t size() const{return window<void>::size()/sizeof(T);}
};
@ -238,16 +231,6 @@ class panel{
// friend window;
};
template<class T /*=void*/>
window<T> communicator::make_window(T* t, mpi3::size_t n){
return window<T>(t, n, *this);
}
template<class T /*=void*/>
window<T> communicator::make_window(){
return make_window<T>((T*)nullptr, 0);
}
template<class T> struct reference;
template<class T>
@ -293,30 +276,6 @@ shm_pointer<T> communicator::allocate_shared(MPI_Aint size) const
}
#endif
template<class T>
void communicator::deallocate_shared(pointer<T>){
// MPI_Free_mem(p.base_ptr(rank()));
}
template<class T>
void communicator::deallocate(pointer<T>&, MPI_Aint){
// p.pimpl_->fence();
// MPI_Free_mem(p.local_ptr());
// MPI_Win_free(&p.pimpl_->impl_);
// delete p.pimpl_;
// p.pimpl_ == nullptr;
}
template<class T>
window<T> communicator::make_window(mpi3::size_t size){
mpi3::info inf;
void* ptr;
window<T> ret;
int s = MPI_Win_allocate(size*sizeof(T), sizeof(T), inf.impl_, this->impl_, &ptr, &ret.impl_);
if(s != MPI_SUCCESS) throw std::runtime_error("cannot window_allocate");
return ret;
}
}}
#ifdef _TEST_BOOST_MPI3_WINDOW
@ -329,21 +288,29 @@ using std::cout;
int mpi3::main(int, char*[], mpi3::communicator world){
std::vector<double> darr(world.rank()?0:100);
mpi3::window<double> w = world.make_window(darr.data(), darr.size());
w.fence();
std::vector<double> darr(world.rank()?0:20);
std::iota(darr.begin(), darr.end(), 0);
mpi3::window<double> win{world, darr.data(), darr.size()};
if(world.rank() == 0){
std::cout << win.size() << std::endl;
assert( win.size() == 20 );
assert( win.base()[13] == 13 );
}else{
assert(win.size() == 0);
assert(win.base() == nullptr );
}
win.fence();
if(world.rank() == 0){
std::vector<double> a = {5., 6.};
w.put(a.begin(), a.end(), 0);
win.put(a.begin(), a.end(), 0);
}
world.barrier();
w.fence();
mpi3::communicator{win.get_group()}.barrier();
win.fence();
std::vector<double> b(2);
w.get(b.begin(), b.end(), 0);
w.fence();
assert( b[0] == 5.);
world.barrier();
win.get(b.begin(), b.end(), 0);
win.fence();
assert( b[0] == 5. and b[1] == 6. );
return 0;
}

642
nexus/bin/nxs-test
View File

@ -35,11 +35,6 @@ for d in dirs_check:
# enforce association between nxs-test executable and Nexus library
sys.path.insert(0,library_dir)
# global data
global_data = dict(
job_ref_table = False,
)
from generic import generic_settings,object_interface,NexusError
generic_settings.raise_error = True
@ -1908,643 +1903,52 @@ def settings_operation():
def machines():
# divert logging function
nlog_divert()
nlabel('imports')
from random import randint
from generic import obj
from machines import Machine,Workstation,InteractiveCluster,Supercomputer
from machines import Job
from machines import get_machine,get_machine_name
nunit('import')
nunit('cpu_count')
nlabel('unimplemented_base_interface')
arg0 = None
arg1 = None
try:
Machine.query_queue(arg0)
nfail()
except:
npass()
#end try
try:
Machine.submit_jobs(arg0)
nfail()
except:
npass()
#end try
try:
Machine.process_job(arg0,arg1)
nfail()
except:
npass()
#end try
try:
Machine.process_job_options(arg0,arg1)
nfail()
except:
npass()
#end try
try:
Machine.write_job(arg0,arg1,file=False)
nfail()
except:
npass()
#end try
try:
Machine.submit_job(arg0,arg1)
nfail()
except:
npass()
#end try
nunit('options')
nunit('job_init')
nlabel('file_setup')
nassert(len(Machine.machines)>0)
for m in Machine.machines:
nassert(isinstance(m,Machine))
exists = m.name in Machine.machines
nassert(exists)
nassert(Machine.exists(m.name))
nassert(Machine.is_unique(m))
m.validate()
#end for
nunit('job_time')
nunit('job_set_id')
nlabel('add')
mtest = Machine.machines.first()
nassert(isinstance(mtest,Machine))
try:
Machine.add(mtest)
nfail()
except:
npass()
#end try
try:
Machine.add('my_machine')
nfail()
except:
npass()
#end try
nunit('job_get_machine')
nunit('job_set_environment')
nlabel('get')
m = Machine.get(mtest.name)
nassert(isinstance(m,Machine))
nassert(id(m)==id(mtest))
try:
Machine.get(m)
nfail()
except:
npass()
#end try
try:
Machine.get('some_nonexistant_machine')
nfail()
except:
npass()
#end try
nunit('job_clone')
nunit('job_serial_clone')
nlabel('intantiation')
# test guards against empty/invalid instantiation
try:
Machine()
nfail()
except:
npass()
#end try
try:
Machine(123)
nfail()
except:
npass()
#end try
# test creation of a new machine
test_name = 'test_machine'
nassert(not Machine.exists(test_name))
m = Machine(name=test_name)
nassert(isinstance(m,Machine))
nassert(Machine.exists(m.name))
nassert(Machine.is_unique(m))
m.validate()
nunit('machine_virtuals')
# test guards against multiple instantiation
try:
Machine(name=test_name)
nfail()
except:
npass()
#end try
nunit('machine_list')
# remove test machine
del Machine.machines.test_machine
nassert(not Machine.exists(test_name))
nunit('machine_add')
nunit('machine_get')
# sort known machines
workstations = obj()
supercomputers = obj()
for machine in Machine.machines:
if isinstance(machine,Workstation):
workstations.append(machine)
elif isinstance(machine,Supercomputer):
supercomputers[machine.name] = machine
else:
Machine.class_error('unknown machine type encountered: {0}'.format(machine.__class__.__name__),'check_process_job_idempotency')
#end if
#end for
nunit('machine_instantiation')
nunit('workstation_init')
nlabel('process_job_idempotency')
def check_process_job_idempotency(nw=10,nwj=10,nsj=10,nij=10):
allow_warn = Machine.allow_warnings
Machine.allow_warnings = False
nunit('workstation_scheduling')
not_idempotent = obj()
nunit('supercomputer_init')
# check workstations
nworkstations = nw
if nworkstations is None:
nworkstations=len(workstations)
#end if
nworkstations = min(nworkstations,len(workstations))
njobs = nwj
for nm in range(nworkstations):
if nworkstations<len(workstations):
machine = workstations.select_random() # select machine at random
else:
machine = workstations[nm]
#end if
cores_min = 1
cores_max = machine.cores
processes_min = 1
processes_max = machine.cores
threads_min = 1
threads_max = machine.cores
job_inputs = []
job_inputs_base = []
for nj in range(njobs): # vary cores
cores = randint(cores_min,cores_max)
threads = randint(threads_min,threads_max)
job_inputs_base.append(obj(cores=cores,threads=threads))
#end for
for nj in range(njobs): # vary processes
processes = randint(processes_min,processes_max)
threads = randint(threads_min,threads_max)
job_inputs_base.append(obj(processes=processes,threads=threads))
#end for
job_inputs.extend(job_inputs_base)
for job_input in job_inputs_base: # run in serial
ji = job_input.copy()
ji.serial = True
job_inputs.append(ji)
#end for
# perform idempotency test
machine_idempotent = True
for job_input in job_inputs:
job = Job(machine=machine.name,**job_input)
job2 = obj.copy(job)
machine.process_job(job2)
machine_idempotent &= job==job2
#end for
if not machine_idempotent:
not_idempotent[machine.name] = machine
#end if
#end for
nunit('supercomputer_scheduling')
# check supercomputers
njobs = nsj
small_node_ceiling = 20
nodes_min = 1
cores_min = 1
threads_min = 1
shared_job_inputs = obj(name='some_job',account='some_account')
for machine in supercomputers:
job_inputs = []
job_inputs_base = []
threads_max = 2*machine.cores_per_node
# sample small number of nodes more heavily
nodes_max = min(small_node_ceiling,machine.nodes)
cores_max = min(small_node_ceiling*machine.cores_per_node,machine.cores)
for nj in range(njobs): # nodes alone
nodes = randint(nodes_min,nodes_max)
threads = randint(threads_min,threads_max)
job_input = obj(nodes=nodes,threads=threads,**shared_job_inputs)
job_inputs_base.append(job_input)
#end for
for nj in range(njobs): # cores alone
cores = randint(cores_min,cores_max)
threads = randint(threads_min,threads_max)
job_input = obj(cores=cores,threads=threads,**shared_job_inputs)
job_inputs_base.append(job_input)
#end for
for nj in range(njobs): # nodes and cores
nodes = randint(nodes_min,nodes_max)
cores = randint(cores_min,cores_max)
threads = randint(threads_min,threads_max)
job_input = obj(nodes=nodes,cores=cores,threads=threads,**shared_job_inputs)
job_inputs_base.append(job_input)
#end for
# sample full node set
nodes_max = machine.nodes
cores_max = machine.cores
for nj in range(njobs): # nodes alone
nodes = randint(nodes_min,nodes_max)
threads = randint(threads_min,threads_max)
job_input = obj(nodes=nodes,threads=threads,**shared_job_inputs)
job_inputs_base.append(job_input)
#end for
for nj in range(njobs): # cores alone
cores = randint(cores_min,cores_max)
threads = randint(threads_min,threads_max)
job_input = obj(cores=cores,threads=threads,**shared_job_inputs)
job_inputs_base.append(job_input)
#end for
for nj in range(njobs): # nodes and cores
nodes = randint(nodes_min,nodes_max)
cores = randint(cores_min,cores_max)
threads = randint(threads_min,threads_max)
job_input = obj(nodes=nodes,cores=cores,threads=threads,**shared_job_inputs)
job_inputs_base.append(job_input)
#end for
job_inputs.extend(job_inputs_base)
# now add serial jobs
for job_input in job_inputs_base:
ji = job_input.copy()
ji.serial = True
job_inputs.append(ji)
#end for
# now add local, serial jobs
for job_input in job_inputs_base:
ji = job_input.copy()
ji.serial = True
ji.local = True
job_inputs.append(ji)
#end for
# perform idempotency test
machine_idempotent = True
for job_input in job_inputs:
job = Job(machine=machine.name,**job_input)
job2 = obj.copy(job)
machine.process_job(job2)
#job_idempotent = job==job2
job_idempotent = object_eq(job,job2)
if not job_idempotent:
d,d1,d2 = object_diff(job,job2,full=True)
change = obj(job_before=obj(d1),job_after=obj(d2))
msg = machine.name+'\n'+str(change)
nfail(msg)
#end if
machine_idempotent &= job_idempotent
#end for
if not machine_idempotent:
not_idempotent[machine.name] = machine
#end if
#end for
nunit('process_job')
if len(not_idempotent)>0:
mlist = ''
for name in sorted(not_idempotent.keys()):
mlist+= '\n '+name
#end for
Machine.class_error('\n\nsome machines failed process_job idempotency test:{0}'.format(mlist))
#end if
Machine.class_log('done checking idempotency')
Machine.allow_warnings = allow_warn
#end def check_process_job_idempotency
nunit('job_run_command')
check_process_job_idempotency()
nunit('write_job')
# test Job.run_command
nlabel('job_run_command')
def parse_job_command(command):
tokens = command.replace(':',' ').split()
launcher = tokens[0]
exe = tokens[-1]
args = []
options = obj()
last_option = None
for t in tokens[1:-1]:
if t.startswith('-'):
options[t] = None
last_option = t
elif last_option is not None:
options[last_option] = t
last_option = None
else:
args.append(t)
#end if
#end for
jc = obj(
launcher = launcher,
executable = exe,
args = args,
options = options,
)
return jc
#end def parse_job_command
def job_commands_equal(c1,c2):
jc1 = parse_job_command(c1)
jc2 = parse_job_command(c2)
return object_eq(jc1,jc2)
#end def job_command_equal
job_run_ref = obj({
('amos' , 'n1' ) : 'srun test.x',
('amos' , 'n1_p1' ) : 'srun test.x',
('amos' , 'n2' ) : 'srun test.x',
('amos' , 'n2_t2' ) : 'srun test.x',
('amos' , 'n2_t2_e' ) : 'srun test.x',
('amos' , 'n2_t2_p2' ) : 'srun test.x',
('bluewaters_xe' , 'n1' ) : 'aprun -n 32 test.x',
('bluewaters_xe' , 'n1_p1' ) : 'aprun -n 1 test.x',
('bluewaters_xe' , 'n2' ) : 'aprun -n 64 test.x',
('bluewaters_xe' , 'n2_t2' ) : 'aprun -d 2 -n 32 test.x',
('bluewaters_xe' , 'n2_t2_e' ) : 'aprun -d 2 -n 32 test.x',
('bluewaters_xe' , 'n2_t2_p2' ) : 'aprun -d 2 -n 4 test.x',
('bluewaters_xk' , 'n1' ) : 'aprun -n 16 test.x',
('bluewaters_xk' , 'n1_p1' ) : 'aprun -n 1 test.x',
('bluewaters_xk' , 'n2' ) : 'aprun -n 32 test.x',
('bluewaters_xk' , 'n2_t2' ) : 'aprun -d 2 -n 16 test.x',
('bluewaters_xk' , 'n2_t2_e' ) : 'aprun -d 2 -n 16 test.x',
('bluewaters_xk' , 'n2_t2_p2' ) : 'aprun -d 2 -n 4 test.x',
('cades' , 'n1' ) : 'mpirun -np 36 test.x',
('cades' , 'n1_p1' ) : 'mpirun -np 1 test.x',
('cades' , 'n2' ) : 'mpirun -np 72 test.x',
('cades' , 'n2_t2' ) : 'mpirun -np 36 --npersocket 9 test.x',
('cades' , 'n2_t2_e' ) : 'mpirun -np 36 --npersocket 9 test.x',
('cades' , 'n2_t2_p2' ) : 'mpirun -np 4 --npersocket 1 test.x',
('cetus' , 'n1' ) : 'runjob --np 16 -p 16 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=1 : test.x',
('cetus' , 'n1_p1' ) : 'runjob --np 1 -p 1 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=1 : test.x',
('cetus' , 'n2' ) : 'runjob --np 32 -p 16 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=1 : test.x',
('cetus' , 'n2_t2' ) : 'runjob --np 16 -p 8 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=2 : test.x',
('cetus' , 'n2_t2_e' ) : 'runjob --np 16 -p 8 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=2 ENV_VAR=1 : test.x',
('cetus' , 'n2_t2_p2' ) : 'runjob --np 4 -p 2 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=2 : test.x',
('chama' , 'n1' ) : 'srun test.x',
('chama' , 'n1_p1' ) : 'srun test.x',
('chama' , 'n2' ) : 'srun test.x',
('chama' , 'n2_t2' ) : 'srun test.x',
('chama' , 'n2_t2_e' ) : 'srun test.x',
('chama' , 'n2_t2_p2' ) : 'srun test.x',
('cooley' , 'n1' ) : 'mpirun -np 12 test.x',
('cooley' , 'n1_p1' ) : 'mpirun -np 1 test.x',
('cooley' , 'n2' ) : 'mpirun -np 24 test.x',
('cooley' , 'n2_t2' ) : 'mpirun -np 12 test.x',
('cooley' , 'n2_t2_e' ) : 'mpirun -np 12 test.x',
('cooley' , 'n2_t2_p2' ) : 'mpirun -np 4 test.x',
('cori' , 'n1' ) : 'srun test.x',
('cori' , 'n1_p1' ) : 'srun test.x',
('cori' , 'n2' ) : 'srun test.x',
('cori' , 'n2_t2' ) : 'srun test.x',
('cori' , 'n2_t2_e' ) : 'srun test.x',
('cori' , 'n2_t2_p2' ) : 'srun test.x',
('edison' , 'n1' ) : 'srun test.x',
('edison' , 'n1_p1' ) : 'srun test.x',
('edison' , 'n2' ) : 'srun test.x',
('edison' , 'n2_t2' ) : 'srun test.x',
('edison' , 'n2_t2_e' ) : 'srun test.x',
('edison' , 'n2_t2_p2' ) : 'srun test.x',
('eos' , 'n1' ) : 'aprun -n 16 test.x',
('eos' , 'n1_p1' ) : 'aprun -n 1 test.x',
('eos' , 'n2' ) : 'aprun -n 32 test.x',
('eos' , 'n2_t2' ) : 'aprun -ss -cc numa_node -d 2 -n 16 test.x',
('eos' , 'n2_t2_e' ) : 'aprun -ss -cc numa_node -d 2 -n 16 test.x',
('eos' , 'n2_t2_p2' ) : 'aprun -ss -cc numa_node -d 2 -n 4 test.x',
('jaguar' , 'n1' ) : 'aprun -n 16 test.x',
('jaguar' , 'n1_p1' ) : 'aprun -n 1 test.x',
('jaguar' , 'n2' ) : 'aprun -n 32 test.x',
('jaguar' , 'n2_t2' ) : 'aprun -d 2 -n 16 test.x',
('jaguar' , 'n2_t2_e' ) : 'aprun -d 2 -n 16 test.x',
('jaguar' , 'n2_t2_p2' ) : 'aprun -d 2 -n 4 test.x',
('komodo' , 'n1' ) : 'mpirun -np 12 test.x',
('komodo' , 'n1_p1' ) : 'mpirun -np 1 test.x',
('komodo' , 'n2' ) : 'mpirun -np 24 test.x',
('komodo' , 'n2_t2' ) : 'mpirun -np 12 test.x',
('komodo' , 'n2_t2_e' ) : 'mpirun -np 12 test.x',
('komodo' , 'n2_t2_p2' ) : 'mpirun -np 4 test.x',
('kraken' , 'n1' ) : 'aprun -n 12 test.x',
('kraken' , 'n1_p1' ) : 'aprun -n 1 test.x',
('kraken' , 'n2' ) : 'aprun -n 24 test.x',
('kraken' , 'n2_t2' ) : 'aprun -d 2 -n 12 test.x',
('kraken' , 'n2_t2_e' ) : 'aprun -d 2 -n 12 test.x',
('kraken' , 'n2_t2_p2' ) : 'aprun -d 2 -n 4 test.x',
('lonestar' , 'n1' ) : 'ibrun -n 12 -o 0 test.x',
('lonestar' , 'n1_p1' ) : 'ibrun -n 1 -o 0 test.x',
('lonestar' , 'n2' ) : 'ibrun -n 24 -o 0 test.x',
('lonestar' , 'n2_t2' ) : 'ibrun -n 12 -o 0 test.x',
('lonestar' , 'n2_t2_e' ) : 'ibrun -n 12 -o 0 test.x',
('lonestar' , 'n2_t2_p2' ) : 'ibrun -n 4 -o 0 test.x',
('matisse' , 'n1' ) : 'mpirun -np 16 test.x',
('matisse' , 'n1_p1' ) : 'mpirun -np 1 test.x',
('matisse' , 'n2' ) : 'mpirun -np 32 test.x',
('matisse' , 'n2_t2' ) : 'mpirun -np 16 test.x',
('matisse' , 'n2_t2_e' ) : 'mpirun -np 16 test.x',
('matisse' , 'n2_t2_p2' ) : 'mpirun -np 4 test.x',
('mira' , 'n1' ) : 'runjob --np 16 -p 16 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=1 : test.x',
('mira' , 'n1_p1' ) : 'runjob --np 1 -p 1 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=1 : test.x',
('mira' , 'n2' ) : 'runjob --np 32 -p 16 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=1 : test.x',
('mira' , 'n2_t2' ) : 'runjob --np 16 -p 8 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=2 : test.x',
('mira' , 'n2_t2_e' ) : 'runjob --np 16 -p 8 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=2 ENV_VAR=1 : test.x',
('mira' , 'n2_t2_p2' ) : 'runjob --np 4 -p 2 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=2 : test.x',
('oic5' , 'n1' ) : 'mpirun -np 32 test.x',
('oic5' , 'n1_p1' ) : 'mpirun -np 1 test.x',
('oic5' , 'n2' ) : 'mpirun -np 64 test.x',
('oic5' , 'n2_t2' ) : 'mpirun -np 32 test.x',
('oic5' , 'n2_t2_e' ) : 'mpirun -np 32 test.x',
('oic5' , 'n2_t2_p2' ) : 'mpirun -np 4 test.x',
('redsky' , 'n1' ) : 'srun test.x',
('redsky' , 'n1_p1' ) : 'srun test.x',
('redsky' , 'n2' ) : 'srun test.x',
('redsky' , 'n2_t2' ) : 'srun test.x',
('redsky' , 'n2_t2_e' ) : 'srun test.x',
('redsky' , 'n2_t2_p2' ) : 'srun test.x',
('serrano' , 'n1' ) : 'srun test.x',
('serrano' , 'n1_p1' ) : 'srun test.x',
('serrano' , 'n2' ) : 'srun test.x',
('serrano' , 'n2_t2' ) : 'srun test.x',
('serrano' , 'n2_t2_e' ) : 'srun test.x',
('serrano' , 'n2_t2_p2' ) : 'srun test.x',
('skybridge' , 'n1' ) : 'srun test.x',
('skybridge' , 'n1_p1' ) : 'srun test.x',
('skybridge' , 'n2' ) : 'srun test.x',
('skybridge' , 'n2_t2' ) : 'srun test.x',
('skybridge' , 'n2_t2_e' ) : 'srun test.x',
('skybridge' , 'n2_t2_p2' ) : 'srun test.x',
('solo' , 'n1' ) : 'srun test.x',
('solo' , 'n1_p1' ) : 'srun test.x',
('solo' , 'n2' ) : 'srun test.x',
('solo' , 'n2_t2' ) : 'srun test.x',
('solo' , 'n2_t2_e' ) : 'srun test.x',
('solo' , 'n2_t2_p2' ) : 'srun test.x',
('stampede2' , 'n1' ) : 'ibrun -n 68 -o 0 test.x',
('stampede2' , 'n1_p1' ) : 'ibrun -n 1 -o 0 test.x',
('stampede2' , 'n2' ) : 'ibrun -n 136 -o 0 test.x',
('stampede2' , 'n2_t2' ) : 'ibrun -n 68 -o 0 test.x',
('stampede2' , 'n2_t2_e' ) : 'ibrun -n 68 -o 0 test.x',
('stampede2' , 'n2_t2_p2' ) : 'ibrun -n 4 -o 0 test.x',
('summit' , 'n1' ) : 'jsrun -a 21 -r 2 -b rs -c 21 -d packed -n 2 -g 0 test.x',
('summit' , 'n1_g6' ) : 'jsrun -a 7 -r 6 -b rs -c 7 -d packed -n 6 -g 1 test.x',
('summit' , 'n2' ) : 'jsrun -a 21 -r 2 -b rs -c 21 -d packed -n 4 -g 0 test.x',
('summit' , 'n2_g6' ) : 'jsrun -a 7 -r 6 -b rs -c 7 -d packed -n 12 -g 1 test.x',
('summit' , 'n2_t2' ) : 'jsrun -a 10 -r 2 -b rs -c 20 -d packed -n 4 -g 0 test.x',
('summit' , 'n2_t2_e' ) : 'jsrun -a 10 -r 2 -b rs -c 20 -d packed -n 4 -g 0 test.x',
('summit' , 'n2_t2_e_g6' ) : 'jsrun -a 3 -r 6 -b rs -c 6 -d packed -n 12 -g 1 test.x',
('summit' , 'n2_t2_g6' ) : 'jsrun -a 3 -r 6 -b rs -c 6 -d packed -n 12 -g 1 test.x',
('supermuc' , 'n1' ) : 'mpiexec -n 28 test.x',
('supermuc' , 'n1_p1' ) : 'mpiexec -n 1 test.x',
('supermuc' , 'n2' ) : 'mpiexec -n 56 test.x',
('supermuc' , 'n2_t2' ) : 'mpiexec -n 28 test.x',
('supermuc' , 'n2_t2_e' ) : 'mpiexec -n 28 test.x',
('supermuc' , 'n2_t2_p2' ) : 'mpiexec -n 4 test.x',
('supermucng' , 'n1' ) : 'mpiexec -n 48 test.x',
('supermucng' , 'n1_p1' ) : 'mpiexec -n 1 test.x',
('supermucng' , 'n2' ) : 'mpiexec -n 96 test.x',
('supermucng' , 'n2_t2' ) : 'mpiexec -n 48 test.x',
('supermucng' , 'n2_t2_e' ) : 'mpiexec -n 48 test.x',
('supermucng' , 'n2_t2_p2' ) : 'mpiexec -n 4 test.x',
('taub' , 'n1' ) : 'mpirun -np 12 test.x',
('taub' , 'n1_p1' ) : 'mpirun -np 1 test.x',
('taub' , 'n2' ) : 'mpirun -np 24 test.x',
('taub' , 'n2_t2' ) : 'mpirun -np 12 test.x',
('taub' , 'n2_t2_e' ) : 'mpirun -np 12 test.x',
('taub' , 'n2_t2_p2' ) : 'mpirun -np 4 test.x',
('theta' , 'n1' ) : 'aprun -e OMP_NUM_THREADS=1 -d 1 -cc depth -j 1 -n 64 -N 64 test.x',
('theta' , 'n1_p1' ) : 'aprun -e OMP_NUM_THREADS=1 -d 1 -cc depth -j 1 -n 1 -N 1 test.x',
('theta' , 'n2' ) : 'aprun -e OMP_NUM_THREADS=1 -d 1 -cc depth -j 1 -n 128 -N 64 test.x',
('theta' , 'n2_t2' ) : 'aprun -e OMP_NUM_THREADS=2 -d 2 -cc depth -j 1 -n 64 -N 32 test.x',
('theta' , 'n2_t2_e' ) : 'aprun -e OMP_NUM_THREADS=2 -d 2 -cc depth -j 1 -n 64 -N 32 test.x',
('theta' , 'n2_t2_p2' ) : 'aprun -e OMP_NUM_THREADS=2 -d 2 -cc depth -j 1 -n 4 -N 2 test.x',
('titan' , 'n1' ) : 'aprun -n 16 test.x',
('titan' , 'n1_p1' ) : 'aprun -n 1 test.x',
('titan' , 'n2' ) : 'aprun -n 32 test.x',
('titan' , 'n2_t2' ) : 'aprun -d 2 -n 16 test.x',
('titan' , 'n2_t2_e' ) : 'aprun -d 2 -n 16 test.x',
('titan' , 'n2_t2_p2' ) : 'aprun -d 2 -n 4 test.x',
('tomcat3' , 'n1' ) : 'mpirun -np 64 test.x',
('tomcat3' , 'n1_p1' ) : 'mpirun -np 1 test.x',
('tomcat3' , 'n2' ) : 'mpirun -np 128 test.x',
('tomcat3' , 'n2_t2' ) : 'mpirun -np 64 test.x',
('tomcat3' , 'n2_t2_e' ) : 'mpirun -np 64 test.x',
('tomcat3' , 'n2_t2_p2' ) : 'mpirun -np 4 test.x',
('uno' , 'n1' ) : 'srun test.x',
('uno' , 'n1_p1' ) : 'srun test.x',
('uno' , 'n2' ) : 'srun test.x',
('uno' , 'n2_t2' ) : 'srun test.x',
('uno' , 'n2_t2_e' ) : 'srun test.x',
('uno' , 'n2_t2_p2' ) : 'srun test.x',
('vesta' , 'n1' ) : 'runjob --np 16 -p 16 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=1 : test.x',
('vesta' , 'n1_p1' ) : 'runjob --np 1 -p 1 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=1 : test.x',
('vesta' , 'n2' ) : 'runjob --np 32 -p 16 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=1 : test.x',
('vesta' , 'n2_t2' ) : 'runjob --np 16 -p 8 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=2 : test.x',
('vesta' , 'n2_t2_e' ) : 'runjob --np 16 -p 8 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=2 ENV_VAR=1 : test.x',
('vesta' , 'n2_t2_p2' ) : 'runjob --np 4 -p 2 $LOCARGS --verbose=INFO --envs OMP_NUM_THREADS=2 : test.x',
})
job_inputs_orig = obj(
n1 = obj(nodes=1),
n1_p1 = obj(nodes=1,processes_per_node=1),
n2 = obj(nodes=2),
n2_t2 = obj(nodes=2,threads=2),
n2_t2_p2 = obj(nodes=2,threads=2,processes_per_node=2),
n2_t2_e = obj(nodes=2,threads=2,env=obj(ENV_VAR=1)),
)
for name in sorted(supercomputers.keys()):
m = supercomputers[name]
if m.requires_account:
acc = 'ABC123'
else:
acc = None
#end if
job_inputs = job_inputs_orig
if name=='summit': # exceptional treatment for summit nodes
job_inputs = job_inputs_orig.copy()
jtypes = list(job_inputs.keys())
for jtype in jtypes:
if 'p' in jtype:
del job_inputs[jtype]
else:
jcpu = job_inputs[jtype]
jcpu.gpus = 0
jgpu = jcpu.copy()
jgpu.gpus = 6
job_inputs[jtype+'_g6'] = jgpu
#end if
#end for
#end if
for jtype in sorted(job_inputs.keys()):
job = Job(app_command = 'test.x',
machine = name,
account = acc,
**job_inputs[jtype]
)
command = job.run_command()
if global_data['job_ref_table']:
sname = "'{0}'".format(name)
stype = "'{0}'".format(jtype)
nlog(" ({0:<16} , {1:<16}) : '{2}',".format(sname,stype,command))
continue
#end if
ref_command = job_run_ref[name,jtype]
if not job_commands_equal(command,ref_command):
nfail('Job.run_command for machine "{0}" does not match the reference\njob inputs: {1}\nreference command: {2}\nincorrect command: {3}'.format(name,job_inputs[jtype],ref_command,command))
#end for
#end for
#end for
nlabel('Job.split_nodes')
for name in sorted(supercomputers.keys()):
m = supercomputers[name]
if m.app_launcher=='srun': # no slurm support yet
continue
#end if
if name=='summit': # no summit support
continue
#end if
if m.requires_account:
acc = 'ABC123'
else:
acc = None
#end if
# make a 4 node job
job = Job(app_command = 'test.x',
machine = name,
account = acc,
nodes = 4,
threads = m.cores_per_node,
)
# split the job into 1 and 3 nodes
job1,job2 = job.split_nodes(1)
# get the split run commands
rc = job.run_command()
rc1 = job1.run_command()
rc2 = job2.run_command()
ns = ' {0} '.format(job.nodes)
ns1 = ' {0} '.format(job1.nodes)
ns2 = ' {0} '.format(job2.nodes)
# verify that node count is in each command
nassert(ns in rc )
nassert(ns1 in rc1)
nassert(ns2 in rc2)
# verify that text on either side of node count
# agrees for original and split commands
rcl ,rcr = rc.split(ns,1)
rc1l,rc1r = rc1.split(ns1,1)
rc2l,rc2r = rc2.split(ns2,1)
rcf = rcl+' '+rcr
rc1f = rc1l+' '+rc1r
rc2f = rc2l+' '+rc2r
nassert(job_commands_equal(rcf,rc1f))
nassert(job_commands_equal(rcf,rc2f))
#end for
# restore logging function
nlog_restore()
nunit_all()
#end def machines

118
nexus/lib/machines.py
View File

@ -57,7 +57,7 @@ from imp import load_source
import re,subprocess
def cpu_count():
def cpu_count():
""" Number of virtual or physical CPUs on this system, i.e.
user/real as output by time(1) when called with an optimally scaling
userspace-only program"""
@ -90,6 +90,11 @@ class Options(DevBase):
#end def __init__
def add(self,**kwargs):
self.transfer_from(kwargs)
#end def add
def read(self,options):
nopts = 0
intext = False
@ -119,11 +124,6 @@ class Options(DevBase):
#end def read
def add(self,**kwargs):
self.transfer_from(kwargs)
#end def add
def write(self):
s = ''
for o in self:
@ -135,7 +135,7 @@ class Options(DevBase):
job_defaults = obj(
job_defaults_assign = obj(
name = 'jobname',
type = None,
directory = None,
@ -159,7 +159,6 @@ job_defaults = obj(
hyperthreads = None,
ppn = None,
gpus = None, # number of gpus per node
compiler = None,
serial = False, # run job serially, no mpi
local = False, # run job locally, no queue submission
days = 0,
@ -183,7 +182,10 @@ job_defaults = obj(
default_cpus_per_task = False, # optionally bypass typical nexus processing for supermucng
ntasks_per_core = None,
cpus_per_task = None,
)
# these are not assigned directly
job_defaults_nonassign = obj(
fake = False,
app = None, # name of/path to application
machine = None,
@ -192,8 +194,11 @@ job_defaults = obj(
app_options = None,
run_options = None,
sub_options = None,
skip_machine = False,
)
job_defaults = obj(job_defaults_assign,job_defaults_nonassign)
class Job(NexusCore):
machine = None #default machine if none is specified in settings
@ -210,12 +215,6 @@ class Job(NexusCore):
job_count = 0
intel_compilers = set(['intel','icc','icpc','ifort'])
gnu_compilers = set(['gnu','gcc','g++','gfortran'])
pgi_compilers = set(['pgi'])
cray_compilers = set(['cray'])
@staticmethod
def restore_default_settings():
Job.machine = None
@ -247,15 +246,16 @@ class Job(NexusCore):
kw.set_optional(**job_defaults)
# extract keywords not assigned
app = kw.delete('app')
machine = kw.delete('machine')
options = kw.delete('options')
app_flags = kw.delete('app_flags')
app_options = kw.delete('app_options')
run_options = kw.delete('run_options')
sub_options = kw.delete('sub_options')
env = kw.delete('env')
fake = kw.delete('fake')
app = kw.delete('app')
machine = kw.delete('machine')
options = kw.delete('options')
app_flags = kw.delete('app_flags')
app_options = kw.delete('app_options')
run_options = kw.delete('run_options')
sub_options = kw.delete('sub_options')
env = kw.delete('env')
fake = kw.delete('fake')
skip_machine = kw.delete('skip_machine')
# assign keywords
self.set(**kw)
@ -303,37 +303,30 @@ class Job(NexusCore):
if self.app_props is None:
self.app_props = []
#end if
if self.compiler is not None:
if self.compiler in self.intel_compilers:
self.compiler = 'intel'
elif self.compiler in self.gnu_compilers:
self.compiler = 'gnu'
elif self.compiler in self.pgi_compilers:
self.compiler = 'pgi'
elif self.compiler in self.cray_compilers:
self.compiler = 'cray'
#end if
#end if
if self.serial:
self.cores = 1
self.nodes = None
#end if
if machine is not None:
if skip_machine:
self.machine = machine
#end if
if env is not None:
self.set_environment(**env)
#end if
#check that the machine exists and have it complete the job info
self.process()
else:
if machine is not None:
self.machine = machine
#end if
if env is not None:
self.set_environment(**env)
#end if
#check that the machine exists and have it complete the job info
self.process()
machine = self.get_machine()
self.batch_mode = machine.in_batch_mode()
machine = self.get_machine()
self.batch_mode = machine.in_batch_mode()
if self.bundled_jobs is not None and not machine.batch_capable:
self.error('running batched/bundled jobs on {0} is either not possible or not yet implemented, sorry.'.format(machine.name))
if self.bundled_jobs is not None and not machine.batch_capable:
self.error('running batched/bundled jobs on {0} is either not possible or not yet implemented, sorry.'.format(machine.name))
#end if
#end if
self.normalize_time()
@ -353,6 +346,7 @@ class Job(NexusCore):
#end def process
# test needed
def process_options(self,machine=None):
if machine is None:
machine = self.get_machine()
@ -361,6 +355,7 @@ class Job(NexusCore):
#end def process_options
# test needed
def initialize(self,sim):
self.set_id()
self.identifier = sim.identifier
@ -420,6 +415,7 @@ class Job(NexusCore):
#end def initialize
# test needed
def renew_app_command(self,sim):
if not self.user_app_command:
self.app_command = sim.app_command()
@ -432,6 +428,7 @@ class Job(NexusCore):
#end def set_id
# remove?
def set_processes(self):
if self.processes is None:
self.error('processes should have been set before now\ncontact the developers and have them fix this','Developer')
@ -462,10 +459,6 @@ class Job(NexusCore):
def divert_out_err(self):
self.identifier += '_divert'
#prefix,ext = self.outfile.split('.',1)
#self.outfile = prefix+'_divert.'+ext
#prefix,ext = self.errfile.split('.',1)
#self.errfile = prefix+'_divert.'+ext
#end def divert_out_err
@ -498,28 +491,22 @@ class Job(NexusCore):
#end if
# remove?
def determine_end_status(self,status):
if not nexus_core.generate_only:
self.successful = False # not really implemented yet
#errfile = os.path.join(self.directory,self.errfile)
#if os.path.exists(errfile):
# fobj = open(errfile,'r')
# errors = fobj.read()
# fobj.close()
# if errors.strip()=='':
# self.successful=True
# #end if
##end if
#end if
#end def determine_end_status
# test needed
def write(self,file=False):
machine = self.get_machine()
return machine.write_job(self,file=file)
#end def write
# test needed
def submit(self):
machine = self.get_machine()
machine.add_job(self)
@ -527,6 +514,7 @@ class Job(NexusCore):
#end def submit
# test needed
def reenter_queue(self):
machine = self.get_machine()
machine.requeue_job(self)
@ -856,6 +844,7 @@ class Machine(NexusCore):
def restore_default_settings(self):
self.account = None
self.user = None
self.local_directory = None
self.app_directory = None
self.app_directories = None
@ -1133,6 +1122,7 @@ class Workstation(Machine):
# test needed
class InteractiveCluster(Workstation):
def __init__(self,*args,**kwargs):
@ -1243,7 +1233,7 @@ class Supercomputer(Machine):
sub_launcher = None,
queue_querier = None,
job_remover = None,
name = None
name = None,
):
if name is None:
if self.name is not None:
@ -1372,6 +1362,8 @@ class Supercomputer(Machine):
USUSP = 'suspended',
SSUSP = 'suspended',
)
elif self.queue_querier=='test_query':
None
else:
self.error('ability to query queue with '+self.queue_querier+' has not yet been implemented')
#end if
@ -1379,11 +1371,13 @@ class Supercomputer(Machine):
#end def __init__
# test needed
def interactive_representation(self,cores):
return InteractiveCluster(self,cores)
#end def interactive_representation
# test needed
def requeue_job(self,job):
if isinstance(job,Job):
jid = job.internal_id
@ -1693,6 +1687,11 @@ class Supercomputer(Machine):
#end if
#end if
#end for
elif self.queue_querier=='test_query': # for testing
# pretend that all jobs have finished
for pid in self.processes.keys():
self.system_queue[pid] = 'complete'
#end for
else:
self.error('ability to query queue with '+self.queue_querier+' has not yet been implemented')
#end if
@ -3281,7 +3280,6 @@ class Tomcat3(Supercomputer):
# workstations
for cores in range(1,128+1):
Workstation('ws'+str(cores),cores,'mpirun'),
Workstation('node'+str(cores),cores,'mpirun'),
#end for
# supercomputers and clusters
# nodes sockets cores ram qslots qlaunch qsubmit qstatus qdelete

5
nexus/lib/numerics.py
View File

@ -83,9 +83,10 @@
import sys
import inspect
import numpy as np
from numpy import array,ndarray,zeros,linspace,pi,exp,sqrt,polyfit,polyval
from numpy import sum,abs,arange,empty,sin,cos,dot,atleast_2d,ogrid
from numpy import ones_like,sign,random,cross,prod
from numpy import ones_like,sign,cross,prod
from numpy.linalg import norm
from generic import obj
from developer import unavailable,warn,error
@ -939,7 +940,7 @@ def equilibration_length(x,tail=.5,plot=False,xlim=None,bounces=2,random=True):
bounce[1] = max(bounce[1],bounce[0])
#print len(x),crossings,crossings[1]-crossings[0]+1
if random:
eqlen = bounce[0]+random.randint(bounce[1]-bounce[0]+1)
eqlen = bounce[0]+np.random.randint(bounce[1]-bounce[0]+1)
else:
eqlen = (bounce[0]+bounce[1])//2
#end if

2
nexus/lib/qmcpack_method_analyzers.py
View File

@ -79,7 +79,7 @@ class MethodAnalyzer(QAanalyzer):
complete &= 'opt' in files
#end if
equil = request.equilibration
nblocks_exclude = 0
nblocks_exclude = -1
if isinstance(equil,int):
nblocks_exclude = equil
elif isinstance(equil,(dict,obj)):

1
nexus/lib/qmcpack_quantity_analyzers.py
View File

@ -134,6 +134,7 @@ class DatAnalyzer(QuantityAnalyzer):
if equilibration!=None and nbe==-1:
self.load_data()
nbe = equilibration_length(self.data[equilibration])
assert nbe>=0, 'Number of equilibration blocks is negative.'
self.method_info.nblocks_exclude = nbe
#end if
#end def __init__

232
nexus/lib/structure.py
View File

@ -26,6 +26,100 @@
#! /usr/bin/env python
"""
The :py:mod:`structure` module provides support for atomic structure I/O,
generation, and manipulation.
List of module contents
-----------------------
Read cif file functions:
* :py:func:`read_cif_celldata`
* :py:func:`read_cif_cell`
* :py:func:`read_cif`
Operations on logical conditions:
* :py:func:`equate`
* :py:func:`negate`
Create a Monkhorst-Pack k-point mesh function
* :py:func:`kmesh`
Tile matrix malipulation functions
* :py:func:`reduce_tilematrix`
* :py:func:`tile_magnetization`
Rotate plane function
* :py:func:`rotate_plane`
Trivial filter function
* :py:func:`trivial_filter`
* :py:class:`MaskFilter`
* :py:func:`optimal_tilematrix`
Base class for :py:class:`Structure` class:
* :py:class:`Sobj`
Base class for :py:class:`DefectStructure`, :py:class:`Crystal`, and :py:class:`Jellium` classes:
* :py:class:`Structure`
SeeK-path functions
* :py:func:`\_getseekpath`
* :py:func:`get_conventional_cell`
* :py:func:`get_primitive_cell`
* :py:func:`get_kpath`
* :py:func:`get_symmetry`
* :py:func:`get_structure_with_bands`
* :py:func:`get_band_tiling`
* :py:func:`get_seekpath_full`
Interpolate structures functions
* :py:func:`interpolate_structures`
Animate structures functions
* :py:func:`structure_animation`
Concrete :py:class:`Structure` classes:
* :py:class:`DefectStructure`
* :py:class:`Crystal`
* :py:class:`Jellium`
Structure generation functions:
* :py:func:`generate_cell`
* :py:func:`generate_structure`
* :py:func:`generate_atom_structure`
* :py:func:`generate_dimer_structure`
* :py:func:`generate_trimer_structure`
* :py:func:`generate_jellium_structure`
* :py:func:`generate_crystal_structure`
* :py:func:`generate_defect_structure`
Read structure functions
* :py:func:`read_structure`
Module contents
---------------
"""
import os
import numpy as np
from copy import deepcopy
@ -1457,21 +1551,149 @@ class Structure(Sobj):
self.skew(d)
#end def stretch
# test needed
def skew(self,skew):
def rotate(self,r,rp=None,passive=False,units="radians",check=True):
"""
Arbitrary rotation of the structure.
Parameters
----------
r : `array_like, float, shape (3,3)` or `array_like, float, shape (3)` or `str`
If a 3x3 matrix, then code executes rotation consistent with this matrix --
it is assumed that the matrix acts on a column-major vector (eg, v'=Rv)
If a three-dimensional array, then the operation of the function depends
on the input type of rp in the following ways:
1. If rp is a scalar, then rp is assumed to be an angle and a rotation
of rp is made about the axis defined by r
2. If rp is a vector, then rp is assumed to be an axis and a rotation is made
such that r aligns with rp
3. If rp is a str, then the rotation is such that r aligns with the
axis given by the str ('x', 'y', 'z', 'a0', 'a1', or 'a2')
If a str then the axis, r, is defined by the input label (e.g. 'x', 'y', 'z', 'a1', 'a2', or 'a3')
and the operation of the function depends on the input type of rp in the following
ways (same as above):
1. If rp is a scalar, then rp is assumed to be an angle and a rotation
of rp is made about the axis defined by r
2. If rp is a vector, then rp is assumed to be an axis and a rotation is made
such that r aligns with rp
3. If rp is a str, then the rotation is such that r aligns with the
axis given by the str ('x', 'y', 'z', 'a0', 'a1', or 'a2')
rp : `array_like, float, shape (3), optional` or `str, optional`
If a 3-dimensional vector is given, then rp is assumed to be an axis and a rotation is made
such that the axis r is aligned with rp.
If a str, then rp is assumed to be an angle and a rotation about the axis defined by r
is made by an angle rp
If a str is given, then rp is assumed to be an axis defined by the given label
(e.g. 'x', 'y', 'z', 'a1', 'a2', or 'a3') and a rotation is made such that the axis r
is aligned with rp.
passive : `bool, optional, default False`
If `True`, perform a passive rotation
If `False`, perform an active rotation
units : `str, optional, default "radians"`
Units of rp, if rp is given as an angle (scalar)
check : `bool, optional, default True`
Perform a check to verify rotation matrix is orthogonal
"""
if rp is not None:
dirmap = dict(x=[1,0,0],y=[0,1,0],z=[0,0,1])
if isinstance(r,str):
if r[0]=='a': # r= 'a0', 'a1', or 'a2'
r = self.axes[int(r[1])]
else: # r= 'x', 'y', or 'z'
r = dirmap[r]
#end if
else:
r = array(r,dtype=float)
if len(r.shape)>1:
self.error('r must be given as a 1-d vector or string, if rp is not None')
#end if
#end if
if isinstance(rp,(int,float)):
if units=="radians" or units=="rad":
theta = float(rp)
else:
theta = float(rp)*np.pi/180.0
c = np.cos(theta)
s = np.sin(theta)
else:
if isinstance(rp,str):
if rp[0]=='a': # rp= 'a0', 'a1', or 'a2'
rp = self.axes[int(rp[1])]
else: # rp= 'x', 'y', or 'z'
rp = dirmap[rp]
#end if
else:
rp = array(rp,dtype=float)
#end if
# go from r,rp to r,theta
c = np.dot(r,rp)/np.linalg.norm(r)/np.linalg.norm(rp)
if abs(c-1)<1e-6:
s = 0.0
r = np.array([1,0,0])
else:
s = np.dot(np.cross(r,rp),np.cross(r,rp))/np.linalg.norm(r)/np.linalg.norm(rp)/np.linalg.norm(np.cross(r,rp))
r = np.cross(r,rp)/np.linalg.norm(np.cross(r,rp))
#end if
# make R from r,theta
R = [[ c+r[0]**2.0*(1.0-c), r[0]*r[1]*(1.0-c)-r[2]*s, r[0]*r[2]*(1.0-c)+r[1]*s],
[r[1]*r[0]*(1.0-c)+r[2]*s, c+r[1]**2.0*(1.0-c), r[1]*r[2]*(1.0-c)-r[0]*s],
[r[2]*r[0]*(1.0-c)-r[1]*s, r[2]*r[1]*(1.0-c)+r[0]*s, c+r[2]**2.0*(1.0-c)]]
else:
R = r
#end if
R = array(R,dtype=float)
if passive:
R = R.T
if check:
if not np.allclose(dot(R,R.T),identity(len(R))):
self.error('the function, rotate, must be given an orthogonal matrix')
#end if
#end if
self.matrix_transform(R)
#end def rotate
# test needed
def matrix_transform(self,A):
"""
Arbitrary transformation matrix (column-major).
Parameters
----------
skew : `array_like, float, shape (3,3)`
Transform the structure using the matrix skew. It is assumed that
skew is in column-major form, i.e., it transforms a vector v as
v' = Tv
"""
A = A.T
axinv = inv(self.axes)
axnew = dot(self.axes,skew)
axnew = dot(self.axes,A)
kaxinv = inv(self.kaxes)
kaxnew = dot(self.kaxes,inv(skew).T)
kaxnew = dot(self.kaxes,inv(A).T)
self.pos = dot(dot(self.pos,axinv),axnew)
self.center = dot(dot(self.center,axinv),axnew)
self.kpoints = dot(dot(self.kpoints,kaxinv),kaxnew)
self.axes = axnew
self.kaxes = kaxnew
if self.folded_structure!=None:
self.folded_structure.skew(skew)
self.folded_structure.matrix_transform(A.T)
#end if
#end def matrix_transform
# test needed
def skew(self,skew):
"""
Arbitrary transformation matrix (row-major).
Parameters
----------
skew : `array_like, float, shape (3,3)`
Transform the structure using the matrix skew. It is assumed that
skew is in row-major form, i.e., it transforms a vector v as
v' = vT
"""
self.matrix_transform(skew.T)
#end def skew

17
nexus/lib/testing.py
View File

@ -277,3 +277,20 @@ def restore_nexus_log():
logging_storage.clear()
assert(len(logging_storage)==0)
#end def restore_nexus_log
# declare test failure
# useful inside try/except blocks
def failed(msg='Test failed.'):
assert False,msg
#end def failed
class TestFailed(Exception):
None
#end class TestFailed
global_data = dict(
job_ref_table = False,
)

1815
nexus/tests/unit/test_machines.py Normal file
View File

File diff suppressed because it is too large Load Diff

2
nexus/tests/unit/test_numerics.py
View File

@ -385,6 +385,8 @@ def test_equilibration_length():
eq = equilibration_length(rs,random=False)
assert(eq==52)
eq = equilibration_length(rs,random=True)
#end def test_equilibration_length

64
nexus/tests/unit/test_structure.py
View File

@ -306,6 +306,70 @@ def test_change_units():
def test_rotate():
import numpy as np
ref = get_reference_structures()
s0 = ref.CuO_prim.copy()
s1 = ref.CuO_prim.copy()
# Test the various parameter choices in the case that rp is given
# Perform rotation taking x-axis to x-axis (original positions should be found)
s1.rotate('x','x')
assert(value_eq(s1.pos,s0.pos))
assert(value_eq(s1.axes,s0.axes))
assert(value_eq(s1.kaxes,s0.kaxes))
# Perform active rotation taking x-coords to y-coords
s1.rotate([1,0,0],[0,1,0])
assert(value_eq(s1.pos,np.array([-s0.pos[:,1],s0.pos[:,0],s0.pos[:,2]]).T))
assert(value_eq(s1.axes,np.array([-s0.axes[:,1],s0.axes[:,0],s0.axes[:,2]]).T))
assert(value_eq(s1.kaxes,np.array([-s0.kaxes[:,1],s0.kaxes[:,0],s0.kaxes[:,2]]).T))
# Perform passive rotation taking x-axis to y-axis (original positions should be found)
s1.rotate('x','y',passive=True)
assert(value_eq(s1.pos,s0.pos))
assert(value_eq(s1.axes,s0.axes))
assert(value_eq(s1.kaxes,s0.kaxes))
# Perform active rotation about z-axis by an angle pi/2
s1.rotate('z',np.pi/2.0)
assert(value_eq(s1.pos,np.array([-s0.pos[:,1],s0.pos[:,0],s0.pos[:,2]]).T))
assert(value_eq(s1.axes,np.array([-s0.axes[:,1],s0.axes[:,0],s0.axes[:,2]]).T))
assert(value_eq(s1.kaxes,np.array([-s0.kaxes[:,1],s0.kaxes[:,0],s0.kaxes[:,2]]).T))
# Perform active rotation taking y-coords to x-coords (original positions should be found)
s1.rotate('y',[1,0,0])
assert(value_eq(s1.pos[-1],s0.pos[-1]))
assert(value_eq(s1.axes[-1],s0.axes[-1]))
assert(value_eq(s1.kaxes[-1],s0.kaxes[-1]))
# Perform active rotation taking a0-coords to a2-coords
s1.rotate(s0.axes[0],s0.axes[2])
assert(value_eq(s1.pos[-1],np.array([-2.15536928,3.46035669,0.86507139])))
assert(value_eq(s1.axes[-1],np.array([-3.91292278,3.02549423,-1.35344154])))
assert(value_eq(s1.kaxes[-1],np.array([-0.90768302,0.83458438,-0.15254555])))
# Perform active rotation taking a2-coords to a0-coords (original positions should be found)
s1.rotate('a2','a0')
assert(value_eq(s1.pos,s0.pos))
assert(value_eq(s1.axes,s0.axes))
assert(value_eq(s1.kaxes,s0.kaxes))
# Test the case where rp is not given
# Perform active rotation taking a2-coords to a0-coords
R = [[0.2570157723433977, 0.6326366344635742,-0.7305571719594085],
[0.4370696746690278, 0.5981289557203555, 0.6717230469572912],
[0.8619240060767753,-0.4919478031900122,-0.12277771249328594]]
s1.rotate(R)
assert(value_eq(s1.pos[-1],np.array([-2.15536928,3.46035669,0.86507139])))
assert(value_eq(s1.axes[-1],np.array([-3.91292278,3.02549423,-1.35344154])))
assert(value_eq(s1.kaxes[-1],np.array([-0.90768302,0.83458438,-0.15254555])))
# A final test which places the structure back into its original form
# Perform passive rotation taking a2-coords to a0-coords (original positions should be found)
s1.rotate([-0.5871158698555267, -0.8034668669004766, -0.09867091342903483],1.7050154439645373,passive=True)
assert(value_eq(s1.pos,s0.pos))
assert(value_eq(s1.axes,s0.axes))
assert(value_eq(s1.kaxes,s0.kaxes))
#end def test_change_units
def test_diagonal_tiling():
ref = get_reference_structures()
diag_tilings = [

2
src/QMCWaveFunctions/CMakeLists.txt
View File

@ -28,7 +28,6 @@ SET(WFBASE_SRCS
SPOSetInfo.cpp
SPOSetInputInfo.cpp
SPOSet.cpp
OptimizableSPOSet.cpp
CompositeSPOSet.cpp
HarmonicOscillator/SHOSet.cpp
HarmonicOscillator/SHOSetBuilder.cpp
@ -56,7 +55,6 @@ SET(JASTROW_SRCS
Jastrow/RPAJastrow.cpp
LatticeGaussianProduct.cpp
LatticeGaussianProductBuilder.cpp
OptimizableSPOBuilder.cpp
Fermion/SPOSetProxy.cpp
Fermion/SPOSetProxyForMSD.cpp
)

14
src/QMCWaveFunctions/Fermion/DiracDeterminantBase.h
View File

@ -43,9 +43,9 @@ public:
SPOVTimer(*TimerManager.createTimer("DiracDeterminantBase::spoval", timer_level_fine)),
SPOVGLTimer(*TimerManager.createTimer("DiracDeterminantBase::spovgl", timer_level_fine))
{
Optimizable = Phi->Optimizable;
Optimizable = Phi->isOptimizable();
is_fermionic = true;
ClassName = "DiracDeterminantBase";
ClassName = "DiracDeterminantBase";
registerTimers();
}
@ -76,26 +76,16 @@ public:
virtual inline void checkInVariables(opt_variables_type& active) override
{
Phi->checkInVariables(active);
Phi->checkInVariables(myVars);
}
virtual inline void checkOutVariables(const opt_variables_type& active) override
{
Phi->checkOutVariables(active);
myVars.clear();
myVars.insertFrom(Phi->myVars);
myVars.getIndex(active);
}
virtual void resetParameters(const opt_variables_type& active) override
{
Phi->resetParameters(active);
for (int i = 0; i < myVars.size(); ++i)
{
int ii = myVars.Index[i];
if (ii >= 0)
myVars[i] = active[ii];
}
}
// To be removed with AoS

2
src/QMCWaveFunctions/Fermion/MultiDiracDeterminant.cpp
View File

@ -392,7 +392,7 @@ MultiDiracDeterminant::MultiDiracDeterminant(SPOSetPtr const& spos, int first)
buildTableGradTimer(*TimerManager.createTimer("MultiDiracDeterminant::buildTableGrad")),
ExtraStuffTimer(*TimerManager.createTimer("MultiDiracDeterminant::ExtraStuff"))
{
(spos->Optimizable == true) ? Optimizable = true : Optimizable = false;
(spos->isOptimizable() == true) ? Optimizable = true : Optimizable = false;
ClassName = "MultiDiracDeterminant";
IsCloned = false;

2
src/QMCWaveFunctions/Fermion/SlaterDetBuilder.cpp
View File

@ -589,7 +589,7 @@ bool SlaterDetBuilder::putDeterminant(xmlNodePtr cur, int spin_group)
targetPsi.setndelay(delay_rank);
#endif
slaterdet_0->add(adet, spin_group);
if (psi->Optimizable)
if (psi->isOptimizable())
slaterdet_0->Optimizable = true;
app_log() << std::endl;

4
src/QMCWaveFunctions/Fermion/SlaterDetOpt.cpp
View File

@ -46,8 +46,8 @@ SlaterDetOpt::SlaterDetOpt(ParticleSet& ptcl, SPOSet* spo_ptr, const int up_or_d
resize(m_nel, m_nmo);
m_nlc = Phi->OrbitalSetSize;
m_nb = Phi->BasisSetSize;
m_nlc = Phi->getOrbitalSetSize();
m_nb = Phi->getBasisSetSize();
// make sure we didn't start with a bad m_nlc
check_index_sanity();

6
src/QMCWaveFunctions/MolecularOrbitals/CuspCorr.cpp
View File

@ -25,7 +25,7 @@ void fillRadFunWithPhiBar(int curOrb_,
{
Psi1 = Phi;
Psi2 = Eta;
int norb = Psi1->OrbitalSetSize;
int norb = Psi1->getOrbitalSetSize();
curOrb = curOrb_;
curCenter = curCenter_;
@ -56,7 +56,7 @@ void fillRadFunWithPhi(int curOrb_,
{
Psi1 = Phi;
Psi2 = Eta;
int norb = Psi1->OrbitalSetSize;
int norb = Psi1->getOrbitalSetSize();
curOrb = curOrb_;
curCenter = curCenter_;
Z = Zion;
@ -177,7 +177,7 @@ RealType execute(int curOrb_,
bool writeout = (thisFile != "NULL");
Psi1 = Phi;
Psi2 = Eta;
int norb = Psi1->OrbitalSetSize;
int norb = Psi1->getOrbitalSetSize();
ValueVector_t X(5);
curOrb = curOrb_;
curCenter = curCenter_;

2
src/QMCWaveFunctions/MolecularOrbitals/CuspCorr.h
View File

@ -104,7 +104,7 @@ public:
{
Psi1 = Phi;
Psi2 = Eta;
int norb = Psi1->OrbitalSetSize;
int norb = Psi1->getOrbitalSetSize();
val1.resize(norb);
grad1.resize(norb);

20
src/QMCWaveFunctions/MolecularOrbitals/LCOrbitalSetOpt.h
View File

@ -274,7 +274,7 @@ public:
this->setBasisSet(bs);
// initialize number of molecular orbitals as zero
this->OrbitalSetSize = 0;
setOrbitalSetSize(0);
}
///////////////////////////////////////////////////////////////////////////////////////////////////
@ -291,7 +291,7 @@ public:
this->setSPOSet(spo);
// initialize number of molecular orbitals as zero
this->OrbitalSetSize = 0;
setOrbitalSetSize(0);
}
///////////////////////////////////////////////////////////////////////////////////////////////////
@ -358,7 +358,7 @@ public:
m_spo_set = spo;
// extract the number of single particle orbitals in the basis set
this->BasisSetSize = m_spo_set->OrbitalSetSize;
this->BasisSetSize = m_spo_set->getOrbitalSetSize();
}
///////////////////////////////////////////////////////////////////////////////////////////////////
@ -380,7 +380,7 @@ public:
retval = new LCOrbitalSetOpt(m_basis_set->makeClone(), m_report_level);
retval->C = C;
retval->setOrbitalSetSize(this->OrbitalSetSize);
retval->setOrbitalSetSize(OrbitalSetSize);
retval->init_LCOrbitalSetOpt(0.0);
retval->m_B = m_B;
@ -429,7 +429,7 @@ public:
void setOrbitalSetSize(int norbs)
{
// record the number of linear combinations (i.e. molecular orbitals)
this->OrbitalSetSize = norbs;
OrbitalSetSize = norbs;
app_log() << "LCOrbitalSetOpt finished setOrbitalSetSize with norbs = " << norbs << std::endl;
}
@ -698,7 +698,7 @@ public:
//app_log() << "this->OrbitalSetSize = " << this->OrbitalSetSize << std::endl;
// check sanity
this->check_input_dim("logdet # of columns", "LCOrbitalSetOpt::evaluate_notranspose", logdet.cols(),
this->OrbitalSetSize);
OrbitalSetSize);
if (logdet.cols() != dlogdet.cols() || logdet.cols() != d2logdet.cols())
throw std::runtime_error("logdet, dlogdet, and d2logdet should have the same number of columns in "
"LCOrbitalSetOpt::evaluate_notranspose");
@ -724,7 +724,7 @@ public:
inline void evaluate(const ParticleSet& P, int iat, ValueVector_t& psi)
{
// check input vector dimension
this->check_input_dim("psi", "LCOrbitalSetOpt::evaluate", psi.size(), this->OrbitalSetSize);
this->check_input_dim("psi", "LCOrbitalSetOpt::evaluate", psi.size(), OrbitalSetSize);
// resize temporary arrays if necessary
if (m_temp_g.size() != BasisSetSize)
@ -757,9 +757,9 @@ public:
inline void evaluate(const ParticleSet& P, int iat, ValueVector_t& psi, GradVector_t& dpsi, ValueVector_t& d2psi)
{
// check sanity
this->check_input_dim("d2psi", "LCOrbitalSetOpt::evaluate", d2psi.size(), this->OrbitalSetSize);
this->check_input_dim("dpsi", "LCOrbitalSetOpt::evaluate", dpsi.size(), this->OrbitalSetSize);
this->check_input_dim("psi", "LCOrbitalSetOpt::evaluate", psi.size(), this->OrbitalSetSize);
this->check_input_dim("d2psi", "LCOrbitalSetOpt::evaluate", d2psi.size(), OrbitalSetSize);
this->check_input_dim("dpsi", "LCOrbitalSetOpt::evaluate", dpsi.size(), OrbitalSetSize);
this->check_input_dim("psi", "LCOrbitalSetOpt::evaluate", psi.size(), OrbitalSetSize);
if (psi.size() != dpsi.size() || psi.size() != d2psi.size())
throw std::runtime_error("psi, dpsi, and d2psi vectors must be the same length in LCOrbitalSetOpt::evaluate");

1
src/QMCWaveFunctions/MolecularOrbitals/LCOrbitalSetWithCorrection.cpp
View File

@ -192,7 +192,6 @@ LCOrbitalSet<BS, false>* LCOrbitalSetWithCorrection<BS, false>::clone2LCOrbitalS
newSPO->C = C;
newSPO->Occ.resize(Occ.size());
newSPO->Occ = Occ;
newSPO->className = "LCOrbitalSet";
return newSPO;
}

31
src/QMCWaveFunctions/OptimizableSPOBuilder.cpp
View File

@ -1,31 +0,0 @@
//////////////////////////////////////////////////////////////////////////////////////
// This file is distributed under the University of Illinois/NCSA Open Source License.
// See LICENSE file in top directory for details.
//
// Copyright (c) 2016 Jeongnim Kim and QMCPACK developers.
//
// File developed by: Ken Esler, kpesler@gmail.com, University of Illinois at Urbana-Champaign
// Jeongnim Kim, jeongnim.kim@gmail.com, University of Illinois at Urbana-Champaign
// Jeremy McMinnis, jmcminis@gmail.com, University of Illinois at Urbana-Champaign
// Jaron T. Krogel, krogeljt@ornl.gov, Oak Ridge National Laboratory
//
// File created by: Ken Esler, kpesler@gmail.com, University of Illinois at Urbana-Champaign
//////////////////////////////////////////////////////////////////////////////////////
#include "QMCWaveFunctions/OptimizableSPOBuilder.h"
#include "OhmmsData/AttributeSet.h"
namespace qmcplusplus
{
OptimizableSPOBuilder::OptimizableSPOBuilder(ParticleSet& p, PtclPoolType& psets, Communicate* comm, xmlNodePtr cur)
: SPOSetBuilder(comm), targetPtcl(&p)
{}
SPOSet* OptimizableSPOBuilder::createSPOSetFromXML(xmlNodePtr cur)
{
OptimizableSPOSet* spo = new OptimizableSPOSet();
return spo;
}
} // namespace qmcplusplus

43
src/QMCWaveFunctions/OptimizableSPOBuilder.h
View File

@ -1,43 +0,0 @@
//////////////////////////////////////////////////////////////////////////////////////
// This file is distributed under the University of Illinois/NCSA Open Source License.
// See LICENSE file in top directory for details.
//
// Copyright (c) 2016 Jeongnim Kim and QMCPACK developers.
//
// File developed by: Ken Esler, kpesler@gmail.com, University of Illinois at Urbana-Champaign
// Jeremy McMinnis, jmcminis@gmail.com, University of Illinois at Urbana-Champaign
// Jaron T. Krogel, krogeljt@ornl.gov, Oak Ridge National Laboratory
// Mark A. Berrill, berrillma@ornl.gov, Oak Ridge National Laboratory
//
// File created by: Jeongnim Kim, jeongnim.kim@gmail.com, University of Illinois at Urbana-Champaign
//////////////////////////////////////////////////////////////////////////////////////
#ifndef OPTIMIZABLE_SPO_BUILDER_H
#define OPTIMIZABLE_SPO_BUILDER_H
#include "QMCWaveFunctions/SPOSetBuilder.h"
#include "QMCWaveFunctions/OptimizableSPOSet.h"
namespace qmcplusplus
{
class OptimizableSPOBuilder : public SPOSetBuilder
{
protected:
typedef std::map<std::string, ParticleSet*> PtclPoolType;
typedef std::map<std::string, SPOSet*> SPOPoolType;
ParticleSet* targetPtcl;
public:
OptimizableSPOBuilder(ParticleSet& p, PtclPoolType& psets, Communicate* comm, xmlNodePtr cur);
/** initialize the Antisymmetric wave function for electrons
*@param cur the current xml node
*/
SPOSet* createSPOSetFromXML(xmlNodePtr cur);
// SPOSet* createSPOSetFromXML(xmlNodePtr cur, SPOPool_t& spo_pool);
};
} // namespace qmcplusplus
#endif

618
src/QMCWaveFunctions/OptimizableSPOSet.cpp
View File

@ -1,618 +0,0 @@
//////////////////////////////////////////////////////////////////////////////////////
// This file is distributed under the University of Illinois/NCSA Open Source License.
// See LICENSE file in top directory for details.
//
// Copyright (c) 2016 Jeongnim Kim and QMCPACK developers.
//
// File developed by: Ken Esler, kpesler@gmail.com, University of Illinois at Urbana-Champaign
// Miguel Morales, moralessilva2@llnl.gov, Lawrence Livermore National Laboratory
// Jeongnim Kim, jeongnim.kim@gmail.com, University of Illinois at Urbana-Champaign
// Jeremy McMinnis, jmcminis@gmail.com, University of Illinois at Urbana-Champaign
// Ye Luo, yeluo@anl.gov, Argonne National Laboratory
// Mark A. Berrill, berrillma@ornl.gov, Oak Ridge National Laboratory
//
// File created by: Ken Esler, kpesler@gmail.com, University of Illinois at Urbana-Champaign
//////////////////////////////////////////////////////////////////////////////////////
#include "Message/Communicate.h"
#include "QMCWaveFunctions/OptimizableSPOSet.h"
#include "Numerics/OhmmsBlas.h"
#include "OhmmsData/AttributeSet.h"
namespace qmcplusplus
{
void OptimizableSPOSet::addParameter(std::string id, int iorb, int basis) {}
bool OptimizableSPOSet::put(xmlNodePtr node, SPOPool_t& spo_pool)
{
std::string gsName, basisName;
bool same_k = false;
bool mapped_k = false;
bool same_orbital = false;
bool lapl = false;
derivScale = 1.0;
OhmmsAttributeSet attrib;
attrib.add(gsName, "gs_sposet");
attrib.add(basisName, "basis_sposet");
attrib.add(same_k, "same_k");
attrib.add(mapped_k, "mapped_k");
attrib.add(same_orbital, "same_orbital");
attrib.add(lapl, "lapl");
attrib.add(N, "size");
attrib.add(derivScale, "scale");
attrib.add(thr, "thr");
attrib.put(node);
if (N == 0)
{
app_error() << "You must specify \"size\" attribute for linearopt sposet.\n";
abort();
}
/////////////////////////////////////
// First, find ground-state SPOSet //
/////////////////////////////////////
if (gsName == "")
{
app_error() << "You must supply \"gs_sposet\". Aborting.\n";
abort();
}
SPOPool_t::iterator iter = spo_pool.find(gsName);
if (iter == spo_pool.end())
{
app_error() << "No sposet named \"" << gsName << "\" found. Abort.\n";
abort();
}
else
{
app_log() << " Found ground-state SPOSet \"" << gsName << "\".\n";
GSOrbitals = iter->second;
}
//////////////////////////////////////
// Now, find basis SPOSet from pool //
//////////////////////////////////////
iter = spo_pool.find(basisName);
if (iter != spo_pool.end())
{
BasisOrbitals = iter->second;
app_log() << " Found basis SPOSet \"" << basisName << "\".\n";
}
if (BasisOrbitals == GSOrbitals)
BasisOrbitals = 0;
if (BasisOrbitals)
M = BasisOrbitals->getOrbitalSetSize();
else
M = GSOrbitals->getOrbitalSetSize() - N;
resize(N, M);
app_log() << " linearopt sposet has " << N << " ground-state orbitals and " << M << " basis orbitals.\n";
//Map K points to orbitals.
ValueMatrix_t allowedOrbs;
allowedOrbs.resize(N, M);
// if (mapped_k)
// {
// std::cerr <<"Not Done"<< std::endl;
// }
// else
if (same_k)
{
int off = BasisOrbitals ? 0 : N;
SPOSet* basOrbs = BasisOrbitals ? BasisOrbitals : GSOrbitals;
for (int igs = 0; igs < N; igs++)
{
PosType k_gs = GSOrbitals->get_k(igs);
for (int ib = 0; ib < M; ib++)
{
PosType k_b = basOrbs->get_k(ib + off);
if (dot(k_gs - k_b, k_gs - k_b) < 1.0e-6)
{
// app_log()<<"Allowing ("<<igs<<","<<ib<<") :";
// for(int ix=0;ix<3;ix++) app_log()<<k_gs[ix]<<" ";
// app_log()<<" : ";
// for(int ix=0;ix<3;ix++) app_log()<<k_b[ix]<<" ";
// app_log()<< std::endl;
allowedOrbs(igs, ib) = 1;
}
else
{
// app_log()<<"Not allowing ("<<igs<<","<<ib<<") :";
// for(int ix=0;ix<3;ix++) app_log()<<k_gs[ix]<<" ";
// app_log()<<" : ";
// for(int ix=0;ix<3;ix++) app_log()<<k_b[ix]<<" ";
// app_log()<< std::endl;
allowedOrbs(igs, ib) = 0;
}
}
}
}
else if (mapped_k)
{
xmlNodePtr xmlCoefs = node->xmlChildrenNode;
while (xmlCoefs != NULL)
{
std::string cname((const char*)xmlCoefs->name);
if (cname == "orbitalmap")
{
std::string type("0");
OhmmsAttributeSet cAttrib;
cAttrib.add(type, "type");
cAttrib.put(xmlCoefs);
if (type != "Array")
{
// app_error() << "Unknown correlation type " + type + " in OptimizableSPOSet." + "Resetting to \"Array\"\n";
xmlNewProp(xmlCoefs, (const xmlChar*)"type", (const xmlChar*)"Array");
}
std::vector<RealType> params;
putContent(params, xmlCoefs);
int indx(0);
if (params.size() == N * M)
for (int igs = 0; igs < N; igs++)
{
for (int ib = 0; ib < M; ib++)
{
allowedOrbs(igs, ib) = params[indx++];
}
}
else
{
app_error() << "Map size is incorrect. parameters given= " << params.size() << " parameters needed= " << M * N
<< std::endl;
}
}
xmlCoefs = xmlCoefs->next;
}
}
else if (same_orbital)
{
if (M < N)
{
app_error() << "Map size is incorrect. N<M for same orbital mapping" << std::endl;
}
allowedOrbs = 0;
for (int igs = 0; igs < N; igs++)
allowedOrbs(igs, igs) = 1;
}
else if (lapl)
{
for (int igs = 0; igs < N; igs++)
for (int ib = 0; ib < (M / N); ib++)
allowedOrbs(igs, igs + N * ib) = 1;
}
else
{
for (int igs = 0; igs < N; igs++)
{
for (int ib = 0; ib < M; ib++)
{
allowedOrbs(igs, ib) = 1;
}
}
}
// Now, look for coefficients element
xmlNodePtr xmlCoefs = node->xmlChildrenNode;
while (xmlCoefs != NULL)
{
std::string cname((const char*)xmlCoefs->name);
if (cname == "coefficients")
{
std::string type("0"), id("0");
int state = -1;
int asize(-1);
OhmmsAttributeSet cAttrib;
cAttrib.add(id, "id");
cAttrib.add(type, "type");
cAttrib.add(asize, "size");
cAttrib.add(state, "state");
cAttrib.put(xmlCoefs);
if (state == -1)
{
app_error() << "You must specify the \"state\" attribute in <linearopt>'s <coefficient>.\n";
}
if (type != "Array")
{
// app_error() << "Unknown correlation type " + type + " in OptimizableSPOSet." + "Resetting to \"Array\"\n";
xmlNewProp(xmlCoefs, (const xmlChar*)"type", (const xmlChar*)"Array");
}
std::vector<RealType> params;
putContent(params, xmlCoefs);
app_log() << "Coefficients for state" << state << ":\n";
// std::cerr << "params.size() = " << params.size() << std::endl;
//If params is missized resize and initialize to zero.
if ((params.size() != asize) && (asize > 0))
{
std::vector<RealType> t_params(params);
params.resize(asize, 0.0);
int ipm(0);
if (t_params.size() > 0)
{
for (int ib = 0; ib < M; ib++)
{
#if defined(QMC_COMPLEX)
if ((allowedOrbs(state, ib).real() > 0) && (std::abs(t_params[ipm]) > thr))
{
params[ib] = t_params[ipm++];
#else
if ((allowedOrbs(state, ib) > 0) && (std::abs(t_params[ipm]) > thr))
{
params[ib] = t_params[ipm++];
#endif
}
else
{
params[ib] = 0;
}
}
}
}
else
{
if (params.size())
asize = params.size();
}
// for (int i=0; i< params.size(); i++) {
ValueMatrix_t& cref(*C);
int indx = 0;
for (int i = 0; i < M; i++)
{
std::stringstream sstr;
#ifndef QMC_COMPLEX
// ParamPointers.push_back(&(C(state,i)));
// ParamIndex.push_back(TinyVector<int,2>(state,i));
sstr << id << "_" << indx;
if (allowedOrbs(state, i))
{
ParamPointers.push_back(&(cref(state, i)));
ParamIndex.push_back(TinyVector<int, 2>(state, i));
if ((indx < asize) && (std::abs(params[indx]) > thr))
cref(state, i) = params[indx];
else
cref(state, i) = 0.0;
myVars.insert(sstr.str(), cref(state, i), true, optimize::SPO_P);
indx++;
}
else
{
cref(state, i) = 0.0;
// myVars.insert(sstr.str(),C(state,i),false,optimize::LINEAR_P);
}
#else
//ParamPointers.push_back(&(C(state,i).real()));
//ParamPointers.push_back(&(C(state,i).imag()));
ParamPointers.push_back(reinterpret_cast<RealType*>(cref[state] + i));
ParamPointers.push_back(reinterpret_cast<RealType*>(cref[state] + i) + 1);
ParamIndex.push_back(TinyVector<int, 2>(state, i));
ParamIndex.push_back(TinyVector<int, 2>(state, i));
sstr << id << "_" << 2 * i + 0;
myVars.insert(sstr.str(), cref(state, i).real(), true, optimize::SPO_P);
sstr << id << "_" << 2 * i + 1;
myVars.insert(sstr.str(), cref(state, i).imag(), true, optimize::SPO_P);
#endif
}
for (int i = 0; i < params.size(); i++)
{
char buf[100];
snprintf(buf, 100, " %12.5f\n", params[i]);
app_log() << buf;
}
}
xmlCoefs = xmlCoefs->next;
}
return true;
}
void OptimizableSPOSet::resetTargetParticleSet(ParticleSet& P)
{
GSOrbitals->resetTargetParticleSet(P);
if (BasisOrbitals)
BasisOrbitals->resetTargetParticleSet(P);
}
void OptimizableSPOSet::setOrbitalSetSize(int norbs) { OrbitalSetSize = norbs; }
void OptimizableSPOSet::checkInVariables(opt_variables_type& active) { active.insertFrom(myVars); }
void OptimizableSPOSet::checkOutVariables(const opt_variables_type& active) { myVars.getIndex(active); }
void OptimizableSPOSet::resetParameters(const opt_variables_type& active)
{
for (int i = 0; i < ParamPointers.size(); i++)
{
int loc=myVars.where(i);
if (loc>=0)
if(std::abs(active[loc])>thr) {
myVars[i]=active[loc];
//FIXME the use of "std::real" needs to be removed after the orbital optimization code is made compatible with complex parameters
*(ParamPointers[i]) = std::real( myVars[i] );
}
}
}
// Obsolete
void OptimizableSPOSet::evaluateDerivatives(ParticleSet& P,
int iat,
const opt_variables_type& active,
ValueMatrix_t& d_phi,
ValueMatrix_t& d_lapl_phi)
{
for (int i = 0; i < d_phi.size(); i++)
for (int j = 0; j < N; j++)
d_phi[i][j] = d_lapl_phi[i][j] = ValueType();
// Evaluate basis states
if (BasisOrbitals)
{
BasisOrbitals->evaluate(P, iat, BasisVal);
std::vector<TinyVector<int, 2>>::iterator iter;
std::vector<TinyVector<int, 2>>& act = ActiveBasis[iat];
for (iter = act.begin(); iter != act.end(); iter++)
{
int elem = (*iter)[0];
int param = (*iter)[1];
int loc = myVars.where(param);
if (loc >= 0)
;
}
}
else
{}
}
void OptimizableSPOSet::evaluate(const ParticleSet& P, int iat, ValueVector_t& psi)
{
GSOrbitals->evaluate(P, iat, GSVal);
if (BasisOrbitals)
{
BasisOrbitals->evaluate(P, iat, BasisVal);
BLAS::gemv_trans(N, M, C->data(), &(BasisVal[0]), &(psi[0]));
}
else
BLAS::gemv_trans(N, M, C->data(), &(GSVal[N]), &(psi[0]));
// for (int i=0; i<N; i++) {
// psi[i] = 0.0;
// for (int j=0; j<M; j++)
// psi[i] += C(i,j)*GSVal[N+j];
// }
for (int i = 0; i < N; i++)
psi[i] += GSVal[i];
}
void OptimizableSPOSet::evaluate(const ParticleSet& P,
int iat,
ValueVector_t& psi,
GradVector_t& dpsi,
ValueVector_t& d2psi)
{
GSOrbitals->evaluate(P, iat, GSVal, GSGrad, GSLapl);
const ValueMatrix_t& cref(*C);
if (BasisOrbitals)
{
BasisOrbitals->evaluate(P, iat, BasisVal, BasisGrad, BasisLapl);
BLAS::gemv_trans(N, M, C->data(), &(BasisVal[0]), &(psi[0]));
BLAS::gemv_trans(N, M, C->data(), &(BasisLapl[0]), &(d2psi[0]));
for (int iorb = 0; iorb < N; iorb++)
{
psi[iorb] += GSVal[iorb];
dpsi[iorb] = GSGrad[iorb];
d2psi[iorb] += GSLapl[iorb];
for (int ibasis = 0; ibasis < M; ibasis++)
for (int dim = 0; dim < OHMMS_DIM; dim++)
dpsi[iorb][dim] += cref(iorb, ibasis) * BasisGrad[ibasis][dim];
}
}
else
{
for (int iorb = 0; iorb < N; iorb++)
{
psi[iorb] = GSVal[iorb];
dpsi[iorb] = GSGrad[iorb];
d2psi[iorb] = GSLapl[iorb];
for (int ibasis = 0; ibasis < M; ibasis++)
{
psi[iorb] += cref(iorb, ibasis) * GSVal[N + ibasis];
for (int dim = 0; dim < OHMMS_DIM; dim++)
dpsi[iorb][dim] += cref(iorb, ibasis) * GSGrad[N + ibasis][dim];
d2psi[iorb] += cref(iorb, ibasis) * GSLapl[N + ibasis];
}
}
// BLAS::gemv_trans (N, M, C.data(), &(GSVal[N]), &(psi[0]));
// BLAS::gemv_trans (N, M, C.data(), &(GSLapl[N]), &(d2psi[0]));
}
// for (int i=0; i<N; i++) {
// psi[i] += GSVal[i];
// d2psi[i] += GSLapl[i];
// }
}
void OptimizableSPOSet::evaluateBasis(const ParticleSet& P,
int first,
int last,
ValueMatrix_t& basis_val,
GradMatrix_t& basis_grad,
ValueMatrix_t& basis_lapl)
{
if (BasisOrbitals)
BasisOrbitals->evaluate_notranspose(P, first, last, basis_val, basis_grad, basis_lapl);
else
{
for (int iat = first; iat < last; iat++)
{
GSOrbitals->evaluate(P, iat, GSVal, GSGrad, GSLapl);
for (int i = 0; i < M; i++)
{
basis_val(iat - first, i) = GSVal[N + i];
basis_grad(iat - first, i) = GSGrad[N + i];
basis_lapl(iat - first, i) = GSLapl[N + i];
}
}
}
}
void OptimizableSPOSet::copyParamsFromMatrix(const opt_variables_type& active,
const Matrix<RealType>& mat,
std::vector<RealType>& destVec)
{
for (int ip = 0; ip < myVars.size(); ip++)
{
int loc = myVars.where(ip);
if (loc >= 0)
{
TinyVector<int, 2> idx = ParamIndex[ip];
destVec[loc] += derivScale * mat(idx[0], idx[1]);
}
}
}
void OptimizableSPOSet::copyParamsFromMatrix(const opt_variables_type& active,
const Matrix<ComplexType>& mat,
std::vector<RealType>& destVec)
{
for (int ip = 0; ip < myVars.size(); ip += 2)
{
int loc = myVars.where(ip);
if (loc >= 0)
{
TinyVector<int, 2> idx = ParamIndex[ip];
assert(ParamIndex[ip + 1][0] == idx[0] && ParamIndex[ip + 1][1] == idx[1]);
destVec[loc] += derivScale * mat(idx[0], idx[1]).real();
loc = myVars.where(ip + 1);
assert(loc >= 0);
destVec[loc] += derivScale * mat(idx[0], idx[1]).imag();
}
}
}
void OptimizableSPOSet::evaluate_notranspose(const ParticleSet& P,
int first,
int last,
ValueMatrix_t& logdet,
GradMatrix_t& dlogdet,
ValueMatrix_t& d2logdet)
{
// std::cerr << "GSValMatrix.size =(" << GSValMatrix.size(0) << ", " << GSValMatrix.size(1) << ")\n";
// std::cerr << "GSGradMatrix.size =(" << GSGradMatrix.size(0) << ", " << GSGradMatrix.size(1) << ")\n";
// std::cerr << "GSLaplMatrix.size =(" << GSLaplMatrix.size(0) << ", " << GSLaplMatrix.size(1) << ")\n";
// std::cerr << "first=" << first << " last=" << last << std::endl;
GSOrbitals->evaluate_notranspose(P, first, last, GSValMatrix, GSGradMatrix, GSLaplMatrix);
if (BasisOrbitals)
{
BasisOrbitals->evaluate_notranspose(P, first, last, BasisValMatrix, BasisGradMatrix, BasisLaplMatrix);
//Note to Ken:
//Use Numerics/MatrixOperators.h
//for C=AB MatrixOperators::product(C,BasisValMatrix,logdet);
//for C=AB^t MatrixOperators::ABt(C,BasisValMatrix,logdet);
BLAS::gemm('T', 'N', N, N, M, 1.0, C->data(), M, BasisValMatrix.data(), M, 0.0, logdet.data(), N);
logdet += GSValMatrix;
BLAS::gemm('T', 'N', N, N, M, 1.0, C->data(), M, BasisLaplMatrix.data(), M, 0.0, d2logdet.data(), N);
d2logdet += GSLaplMatrix;
// Gradient part.
for (int dim = 0; dim < OHMMS_DIM; dim++)
{
for (int i = 0; i < M; i++)
for (int j = 0; j < N; j++)
GradTmpSrc(i, j) = BasisGradMatrix(i, j)[dim];
BLAS::gemm('T', 'N', N, N, M, 1.0, C->data(), M, GradTmpSrc.data(), M, 0.0, GradTmpDest.data(), N);
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
dlogdet(i, j)[dim] = GradTmpDest(i, j) + GSGradMatrix(i, j)[dim];
}
}
else
{
// HACK HACK HACK
// BLAS::gemm ('T', 'N', N, N, M, 1.0, C.data(),
// M, GSValMatrix.data()+N, M+N, 0.0, logdet.data(), N);
// for (int i=0; i<N; i++)
// for (int j=0; j<N; j++)
// logdet(i,j) += GSValMatrix(i,j);
for (int iel = 0; iel < N; iel++)
for (int iorb = 0; iorb < N; iorb++)
{
logdet(iel, iorb) = GSValMatrix(iel, iorb);
const ValueType* restrict cptr = (*C)[iorb];
for (int ibasis = 0; ibasis < M; ibasis++)
logdet(iel, iorb) += cptr[ibasis] * GSValMatrix(iel, N + ibasis);
}
// logdet += GSValMatrix;
// BLAS::gemm ('T', 'N', N, N, M, 1.0, C.data(),
// M, GSLaplMatrix.data()+N, M+N, 0.0, d2logdet.data(), N);
for (int iel = 0; iel < N; iel++)
for (int iorb = 0; iorb < N; iorb++)
{
const ValueType* restrict cptr = (*C)[iorb];
d2logdet(iel, iorb) = GSLaplMatrix(iel, iorb);
for (int ibasis = 0; ibasis < M; ibasis++)
d2logdet(iel, iorb) += cptr[ibasis] * GSLaplMatrix(iel, N + ibasis);
}
//d2logdet += GSLaplMatrix;
// Gradient part.
for (int i = 0; i < N; i++)
for (int iorb = 0; iorb < N; iorb++)
{
const ValueType* restrict cptr = (*C)[iorb];
dlogdet(i, iorb) = GSGradMatrix(i, iorb);
for (int n = 0; n < M; n++)
dlogdet(i, iorb) += cptr[n] * GSGradMatrix(i, N + n);
}
// for (int dim=0; dim<OHMMS_DIM; dim++) {
// for (int i=0; i<M; i++)
// for (int j=0; j<N; j++)
// GradTmpSrc(i,j) = GSGradMatrix(i,j+N)[dim];
// BLAS::gemm ('T', 'N', N, N, M, 1.0, C.data(), M,
// GradTmpSrc.data(), M, 0.0, GradTmpDest.data(), N);
// for (int i=0; i<N; i++)
// for (int j=0; j<N; j++)
// dlogdet(i,j)[dim] = GradTmpDest(i,j) + GSGradMatrix(i,j)[dim];
//}
}
}
SPOSet* OptimizableSPOSet::makeClone() const
{
SPOSet *gs, *basis(0);
OptimizableSPOSet* clone;
gs = GSOrbitals->makeClone();
if (BasisOrbitals)
basis = BasisOrbitals->makeClone();
clone = new OptimizableSPOSet(N, gs, basis);
clone->IsCloned = true;
clone->C = C; //just pointer
clone->thr = thr;
clone->myVars = myVars;
clone->derivScale = derivScale;
#if 0
clone->ParamPointers.clear();
clone->ParamIndex=ParamIndex;
for(int i=0; i<ParamIndex.size() ; i++)
{
#ifndef QMC_COMPLEX
clone->ParamPointers.push_back(&(clone->C(ParamIndex[i][0],ParamIndex[i][1])));
#else
int ci=ParamIndex[i][0];
int cj=ParamIndex[i][1];
//clone->ParamPointers.push_back(&(clone->C(ParamIndex[i][0],ParamIndex[i][1]).real()));
//clone->ParamPointers.push_back(&(clone->C(ParamIndex[i][0],ParamIndex[i][1]).imag()));
clone->ParamPointers.push_back(reinterpret_cast<RealType*>(clone->C[ci]+cj));
clone->ParamPointers.push_back(reinterpret_cast<RealType*>(clone->C[ci]+cj)+1);
#endif
}
#endif
// for (int i=0; i< N; i++) {
// for (int j=0; j< M; j++) {
// std::stringstream sstr;
// #ifndef QMC_COMPLEX
// clone->ParamPointers.push_back(&(clone->C(i,j)));
// clone->ParamIndex.push_back(TinyVector<int,2>(i,j));
// #else
// clone->ParamPointers.push_back(&(clone->C(i,j).real()));
// clone->ParamPointers.push_back(&(clone->C(i,j).imag()));
// clone->ParamIndex.push_back(TinyVector<int,2>(i,j));
// clone->ParamIndex.push_back(TinyVector<int,2>(i,j));
// #endif
// }
// }
return clone;
}
} // namespace qmcplusplus

226
src/QMCWaveFunctions/OptimizableSPOSet.h
View File

@ -1,226 +0,0 @@
//////////////////////////////////////////////////////////////////////////////////////
// This file is distributed under the University of Illinois/NCSA Open Source License.
// See LICENSE file in top directory for details.
//
// Copyright (c) 2016 Jeongnim Kim and QMCPACK developers.
//
// File developed by: Jeongnim Kim, jeongnim.kim@gmail.com, University of Illinois at Urbana-Champaign
// Ken Esler, kpesler@gmail.com, University of Illinois at Urbana-Champaign
// Miguel Morales, moralessilva2@llnl.gov, Lawrence Livermore National Laboratory
// Jeremy McMinnis, jmcminis@gmail.com, University of Illinois at Urbana-Champaign
// Mark A. Berrill, berrillma@ornl.gov, Oak Ridge National Laboratory
//
// File created by: Jeongnim Kim, jeongnim.kim@gmail.com, University of Illinois at Urbana-Champaign
//////////////////////////////////////////////////////////////////////////////////////
#ifndef OPTIMIZABLE_SPO_SET_H
#define OPTIMIZABLE_SPO_SET_H
#include "QMCWaveFunctions/SPOSet.h"
#include "Numerics/OptimizableFunctorBase.h"
#include "Optimize/VariableSet.h"
namespace qmcplusplus
{
class OptimizableSPOSet : public SPOSet
{
protected:
typedef optimize::VariableSet opt_variables_type;
///typedef for name-value lists
typedef optimize::VariableSet::variable_map_type variable_map_type;
// Number of occupied states, number of basis states
int N, M;
///number of basis
IndexType BasisSetSize;
/** pointer to matrix containing the coefficients
*
* makeClone makes a shallow copy
*/
ValueMatrix_t* C;
///if true, do not clean up
bool IsCloned;
RealType derivScale;
RealType thr;
// If BasisOrbitals==NULL, only GSOrbitals is used and it's evaluate
// functions should return N+M orbitals.
SPOSet *GSOrbitals, *BasisOrbitals;
// The first index is the orbital to be optimized, the second is the
// basis element
// ValueMatrix_t OptCoefs;
// This maps the parameter indices to elements in the C
// matrix. Using pointers allows the use of complex C
// while maintaining real optimizable parameters.
std::vector<RealType*> ParamPointers;
// Maps the parameter index in myVars to an index in the C array
std::vector<TinyVector<int, 2>> ParamIndex;
void addParameter(std::string id, int iorb, int basis);
ValueVector_t GSVal, BasisVal, GSLapl, BasisLapl;
GradVector_t GSGrad, BasisGrad;
ValueMatrix_t GSValMatrix, BasisValMatrix, GSLaplMatrix, BasisLaplMatrix, GradTmpSrc, GradTmpDest;
GradMatrix_t GSGradMatrix, BasisGradMatrix;
// Cache the positions to avoid recomputing GSVal, BasisVal, etc.
// if unnecessary.
std::vector<PosType> CachedPos;
public:
///set of variables to be optimized; These are mapped to the
///C matrix. Update: Moved to SPOSet
// opt_variables_type myVars;
// For each occupied orbital, this lists which of the M
// basis functions are used to construct the optimal orbital.
// This must be initialized by the SPOSet builder. If it is not
// initialized, we assume that C is fully dense.
// The first element of the TinyVector is the basis element
// number. The second is the corresponding parameter number.
std::vector<std::vector<TinyVector<int, 2>>> ActiveBasis;
OptimizableSPOSet() : N(0), M(0), derivScale(10.0), thr(0.0), GSOrbitals(0), BasisOrbitals(0)
{
Optimizable = true;
className = "OptimizableSPOSet";
}
OptimizableSPOSet(int num_orbs, SPOSet* gsOrbs, SPOSet* basisOrbs = 0)
: GSOrbitals(gsOrbs), BasisOrbitals(basisOrbs), derivScale(10.0)
{
N = num_orbs;
setOrbitalSetSize(N);
if (BasisOrbitals)
M = BasisOrbitals->getOrbitalSetSize();
else
M = GSOrbitals->getOrbitalSetSize() - N;
resize(N, M);
Optimizable = true;
}
// bool put(xmlNodePtr cur, SPOPool_t &spo_pool);
void resetTargetParticleSet(ParticleSet& P);
void setOrbitalSetSize(int norbs);
void set_active_basis(std::vector<std::vector<int>>& active)
{
//ActiveBasis = active;
}
// Read coefficients, or create the XML element.
bool put(xmlNodePtr node, SPOPool_t& spo_pool);
// This stores new orbital coefficients int C
void resetParameters(const opt_variables_type& optvars);
// Evaluate the derivative of the optimized orbitals with
// respect to the parameters
void evaluateDerivatives(ParticleSet& P,
int iat,
const opt_variables_type& active,
ValueMatrix_t& d_phi,
ValueMatrix_t& d_lapl_phi);
void checkInVariables(opt_variables_type& active);
void checkOutVariables(const opt_variables_type& active);
// Evaluate routines. These call GSOrbitals->evaluate and possibly
// BasisOrbitals->evaluate, then does the matrix product with
// C.
void evaluate(const ParticleSet& P, int iat, ValueVector_t& psi);
void evaluate(const ParticleSet& P, int iat, ValueVector_t& psi, GradVector_t& dpsi, ValueVector_t& d2psi);
void evaluate_notranspose(const ParticleSet& P,
int first,
int last,
ValueMatrix_t& logdet,
GradMatrix_t& dlogdet,
ValueMatrix_t& d2logdet);
void evaluate_notranspose(const ParticleSet& P,
int first,
int last,
ValueMatrix_t& logdet,
GradMatrix_t& dlogdet,
HessMatrix_t& grad_grad_logdet)
{
APP_ABORT("Need specialization of OptimizableOrbitalSet::evaluate_notranspose() for grad_grad_logdet. \n");
}
void evaluate_notranspose(const ParticleSet& P,
int first,
int last,
ValueMatrix_t& logdet,
GradMatrix_t& dlogdet,
HessMatrix_t& grad_grad_logdet,
GGGMatrix_t& grad_grad_grad_logdet)
{
APP_ABORT("Need specialization of OptimazableOrbitalSet::evaluate_notranspose() for grad_grad_grad_logdet. \n");
}
void evaluateBasis(const ParticleSet& P,
int first,
int last,
ValueMatrix_t& basis_val,
GradMatrix_t& basis_grad,
ValueMatrix_t& basis_lapl);
void copyParamsFromMatrix(const opt_variables_type& active,
const Matrix<RealType>& mat,
std::vector<RealType>& destVec);
void copyParamsFromMatrix(const opt_variables_type& active,
const Matrix<ComplexType>& mat,
std::vector<RealType>& destVec);
void resize(int n, int m)
{
N = n;
M = m;
if (BasisOrbitals)
{
GSVal.resize(N);
GSGrad.resize(N);
GSLapl.resize(N);
BasisVal.resize(M);
BasisGrad.resize(M);
BasisLapl.resize(M);
GSValMatrix.resize(N, N);
GSGradMatrix.resize(N, N);
GSLaplMatrix.resize(N, N);
BasisValMatrix.resize(M, N);
BasisGradMatrix.resize(M, N);
BasisLaplMatrix.resize(M, N);
}
else
{
GSVal.resize(N + M);
GSGrad.resize(N + M);
GSLapl.resize(N + M);
GSValMatrix.resize(N, N + M);
GSGradMatrix.resize(N, N + M);
GSLaplMatrix.resize(N, N + M);
}
GradTmpSrc.resize(M, N);
GradTmpDest.resize(N, N);
if (C == nullptr)
C = new ValueMatrix_t(N, M);
else
C->resize(N, M);
ActiveBasis.resize(N);
BasisSetSize = M;
}
// Make a copy of myself
SPOSet* makeClone() const;
};
} // namespace qmcplusplus
#endif

7
src/QMCWaveFunctions/SPOSet.cpp
View File

@ -26,11 +26,10 @@
namespace qmcplusplus
{
SPOSet::SPOSet()
SPOSet::SPOSet(bool ion_deriv, bool optimizable)
: OrbitalSetSize(0),
Optimizable(false),
ionDerivs(false),
builder_index(-1)
ionDerivs(ion_deriv),
Optimizable(optimizable)
#if !defined(ENABLE_SOA)
,
Identity(false),

27
src/QMCWaveFunctions/SPOSet.h
View File

@ -62,16 +62,6 @@ public:
typedef ParticleSet::Walker_t Walker_t;
typedef std::map<std::string, SPOSet*> SPOPool_t;
///index in the builder list of sposets
int builder_index;
///true if SPO is optimizable
bool Optimizable;
///number of Single-particle orbitals
IndexType OrbitalSetSize;
/// Optimizable variables
opt_variables_type myVars;
///name of the class
std::string className;
/** name of the object
*
* Several user classes can own SPOSet and use objectName as counter
@ -96,7 +86,7 @@ public:
#endif
/** constructor */
SPOSet();
SPOSet(bool ion_deriv = false, bool optimizable = false);
/** destructor
*
@ -110,6 +100,9 @@ public:
#endif
}
// accessor function to Optimizable
inline bool isOptimizable() const { return Optimizable; }
/** return the size of the orbital set
* Ye: this needs to be replaced by getOrbitalSetSize();
*/
@ -472,8 +465,18 @@ protected:
bool putFromXML(xmlNodePtr coeff_ptr);
bool putFromH5(const std::string& fname, xmlNodePtr coeff_ptr);
#endif
protected:
///true, if the derived class has non-zero ionic derivatives.
bool ionDerivs;
const bool ionDerivs;
///true if SPO is optimizable
const bool Optimizable;
///number of Single-particle orbitals
IndexType OrbitalSetSize;
/// Optimizable variables
opt_variables_type myVars;
///name of the class
std::string className;
};
#if defined(ENABLE_SMARTPOINTER)

1
src/QMCWaveFunctions/SPOSetBuilder.cpp
View File

@ -53,7 +53,6 @@ SPOSet* SPOSetBuilder::createSPOSet(xmlNodePtr cur)
if (sposet)
{
//sposet->put(cur); //initialize C and other internal containers
sposet->builder_index = sposets.size();
sposets.push_back(sposet);
}
else

6
src/QMCWaveFunctions/SPOSetBuilderFactory.cpp
View File

@ -36,7 +36,6 @@
#endif
#endif
#include "QMCWaveFunctions/CompositeSPOSet.h"
#include "QMCWaveFunctions/OptimizableSPOBuilder.h"
#include "Utilities/ProgressReportEngine.h"
#include "Utilities/IteratorUtility.h"
#include "OhmmsData/AttributeSet.h"
@ -181,11 +180,6 @@ SPOSetBuilder* SPOSetBuilderFactory::createSPOSetBuilder(xmlNodePtr rootNode)
app_log() << "Harmonic Oscillator SPO set" << std::endl;
bb = new SHOSetBuilder(targetPtcl, myComm);
}
else if (type == "linearopt")
{
//app_log()<<"Optimizable SPO set"<< std::endl;
bb = new OptimizableSPOBuilder(targetPtcl, ptclPool, myComm, rootNode);
}
#if OHMMS_DIM == 3
else if (type.find("spline") < type.size())
{

8
src/QMCWaveFunctions/lcao/CuspCorrection.cpp
View File

@ -163,8 +163,8 @@ void splitPhiEta(int center, const std::vector<bool>& corrCenter, LCAOrbitalSet&
Phi.myBasisSet->queryOrbitalsForSType(correct_this_center, is_s_orbital);
int nOrbs = Phi.OrbitalSetSize;
int bss = Phi.BasisSetSize;
int nOrbs = Phi.getOrbitalSetSize();
int bss = Phi.getBasisSetSize();
for (int i = 0; i < bss; i++)
{
@ -191,8 +191,8 @@ void removeSTypeOrbitals(const std::vector<bool>& corrCenter, LCAOrbitalSet& Phi
Phi.myBasisSet->queryOrbitalsForSType(corrCenter, is_s_orbital);
int nOrbs = Phi.OrbitalSetSize;
int bss = Phi.BasisSetSize;
int nOrbs = Phi.getOrbitalSetSize();
int bss = Phi.getBasisSetSize();
for (int i = 0; i < bss; i++)
{

2
src/QMCWaveFunctions/lcao/CuspCorrection.h
View File

@ -114,7 +114,7 @@ public:
: targetPtcl(targetP), sourcePtcl(sourceP), curOrb(0), curCenter(0)
{
Psi1 = Phi;
int norb = Psi1->OrbitalSetSize;
int norb = Psi1->getOrbitalSetSize();
val1.resize(norb);
grad1.resize(norb);
lap1.resize(norb);

24
src/QMCWaveFunctions/lcao/CuspCorrectionConstruction.cpp
View File

@ -33,14 +33,14 @@ void applyCuspCorrection(const Matrix<CuspCorrectionParameters>& info,
ScopedTimer cuspApplyTimerWrapper(cuspApplyTimer);
LCAOrbitalSet phi = LCAOrbitalSet(lcwc.myBasisSet);
phi.setOrbitalSetSize(lcwc.OrbitalSetSize);
phi.BasisSetSize = lcwc.BasisSetSize;
LCAOrbitalSet phi = LCAOrbitalSet(lcwc.myBasisSet, lcwc.isOptimizable());
phi.setOrbitalSetSize(lcwc.getOrbitalSetSize());
phi.BasisSetSize = lcwc.getBasisSetSize();
phi.setIdentity(false);
LCAOrbitalSet eta = LCAOrbitalSet(lcwc.myBasisSet);
eta.setOrbitalSetSize(lcwc.OrbitalSetSize);
eta.BasisSetSize = lcwc.BasisSetSize;
LCAOrbitalSet eta = LCAOrbitalSet(lcwc.myBasisSet, lcwc.isOptimizable());
eta.setOrbitalSetSize(lcwc.getOrbitalSetSize());
eta.BasisSetSize = lcwc.getBasisSetSize();
eta.setIdentity(false);
@ -201,14 +201,14 @@ void generateCuspInfo(int orbital_set_size,
ScopedTimer createCuspTimerWrapper(cuspCreateTimer);
LCAOrbitalSet phi = LCAOrbitalSet(lcwc.myBasisSet);
phi.setOrbitalSetSize(lcwc.OrbitalSetSize);
phi.BasisSetSize = lcwc.BasisSetSize;
LCAOrbitalSet phi = LCAOrbitalSet(lcwc.myBasisSet, lcwc.isOptimizable());
phi.setOrbitalSetSize(lcwc.getOrbitalSetSize());
phi.BasisSetSize = lcwc.getBasisSetSize();
phi.setIdentity(false);
LCAOrbitalSet eta = LCAOrbitalSet(lcwc.myBasisSet);
eta.setOrbitalSetSize(lcwc.OrbitalSetSize);
eta.BasisSetSize = lcwc.BasisSetSize;
LCAOrbitalSet eta = LCAOrbitalSet(lcwc.myBasisSet, lcwc.isOptimizable());
eta.setOrbitalSetSize(lcwc.getOrbitalSetSize());
eta.BasisSetSize = lcwc.getBasisSetSize();
eta.setIdentity(false);

18
src/QMCWaveFunctions/lcao/LCAOrbitalBuilder.cpp
View File

@ -443,15 +443,19 @@ SPOSet* LCAOrbitalBuilder::createSPOSetFromXML(xmlNodePtr cur)
if (myBasisSet == nullptr)
PRE.error("Missing basisset.", true);
if (optimize == "yes")
app_log() << " SPOSet " << spo_name << " is optimizable\n";
LCAOrbitalSet* lcos = nullptr;
#if !defined(QMC_COMPLEX)
LCAOrbitalSetWithCorrection* lcwc = nullptr;
if (doCuspCorrection)
lcos = lcwc = new LCAOrbitalSetWithCorrection(sourcePtcl, targetPtcl, myBasisSet);
lcos = lcwc = new LCAOrbitalSetWithCorrection(sourcePtcl, targetPtcl, myBasisSet, optimize == "yes");
else
lcos = new LCAOrbitalSet(myBasisSet);
lcos = new LCAOrbitalSet(myBasisSet, optimize == "yes");
#else
lcos = new LCAOrbitalSet(myBasisSet);
lcos = new LCAOrbitalSet(myBasisSet, optimize == "yes");
#endif
loadMO(*lcos, cur);
@ -464,7 +468,7 @@ SPOSet* LCAOrbitalBuilder::createSPOSetFromXML(xmlNodePtr cur)
if (id == "")
id = spo_name;
int orbital_set_size = lcos->OrbitalSetSize;
int orbital_set_size = lcos->getOrbitalSetSize();
Matrix<CuspCorrectionParameters> info(num_centers, orbital_set_size);
bool valid = false;
@ -494,12 +498,6 @@ SPOSet* LCAOrbitalBuilder::createSPOSetFromXML(xmlNodePtr cur)
}
#endif
if (optimize == "yes")
{
lcos->Optimizable = true;
app_log() << " SPOSet " << spo_name << " is optimizable\n";
}
return lcos;
}

6
src/QMCWaveFunctions/lcao/LCAOrbitalSet.cpp
View File

@ -16,11 +16,9 @@
namespace qmcplusplus
{
LCAOrbitalSet::LCAOrbitalSet(basis_type* bs)
: myBasisSet(nullptr), params_supplied(false), BasisSetSize(0), Identity(true)
LCAOrbitalSet::LCAOrbitalSet(basis_type* bs, bool optimize)
: SPOSet(true, optimize), myBasisSet(nullptr), params_supplied(false), BasisSetSize(0), Identity(true)
{
//This SPOSet has an explicit ion dependence, so set this flag.
ionDerivs=true;
if (bs != nullptr)
setBasisSet(bs);
}

8
src/QMCWaveFunctions/lcao/LCAOrbitalSet.h
View File

@ -73,7 +73,7 @@ public:
/** constructor
* @param bs pointer to the BasisSet
*/
LCAOrbitalSet(basis_type* bs = nullptr);
LCAOrbitalSet(basis_type* bs, bool optimize);
LCAOrbitalSet(const LCAOrbitalSet& in) = default;
@ -109,16 +109,11 @@ public:
const size_t NP2,
const std::vector<std::vector<int>>& lookup_tbl);
void checkInVariables(opt_variables_type& active)
{
if (Optimizable)
{
if (myVars.size())
active.insertFrom(myVars);
else
Optimizable = false;
}
}
void checkOutVariables(const opt_variables_type& active)
@ -127,7 +122,6 @@ public:
myVars.getIndex(active);
}
///reset
void resetParameters(const opt_variables_type& active)
{

4
src/QMCWaveFunctions/lcao/LCAOrbitalSetWithCorrection.cpp
View File

@ -15,8 +15,8 @@
namespace qmcplusplus
{
LCAOrbitalSetWithCorrection::LCAOrbitalSetWithCorrection(ParticleSet& ions, ParticleSet& els, basis_type* bs)
: LCAOrbitalSet(bs), cusp(ions, els)
LCAOrbitalSetWithCorrection::LCAOrbitalSetWithCorrection(ParticleSet& ions, ParticleSet& els, basis_type* bs, bool optimize)
: LCAOrbitalSet(bs, optimize), cusp(ions, els)
{}
void LCAOrbitalSetWithCorrection::setOrbitalSetSize(int norbs)

2
src/QMCWaveFunctions/lcao/LCAOrbitalSetWithCorrection.h
View File

@ -34,7 +34,7 @@ struct LCAOrbitalSetWithCorrection : public LCAOrbitalSet
* @param bs pointer to the BasisSet
* @param rl report level
*/
LCAOrbitalSetWithCorrection(ParticleSet& ions, ParticleSet& els, basis_type* bs = nullptr);
LCAOrbitalSetWithCorrection(ParticleSet& ions, ParticleSet& els, basis_type* bs, bool optimize);
LCAOrbitalSetWithCorrection(const LCAOrbitalSetWithCorrection& in) = default;

4
src/QMCWaveFunctions/tests/test_cusp_corr.cpp
View File

@ -108,14 +108,14 @@ TEST_CASE("CuspCorrection He", "[wavefunction]")
typedef NGOBuilder::CenteredOrbitalType COT;
OrbType bs_phi(lcob->myBasisSet);
bs_phi.setOrbitalSetSize(lcob->OrbitalSetSize);
bs_phi.setOrbitalSetSize(lcob->getOrbitalSetSize());
bs_phi.BasisSetSize = lcob->BasisSetSize;
bs_phi.setIdentity(false);
*(bs_phi.C) = *(lcob->C);
OrbType bs_eta(lcob->myBasisSet);
bs_eta.setOrbitalSetSize(lcob->OrbitalSetSize);
bs_eta.setOrbitalSetSize(lcob->getOrbitalSetSize());
bs_eta.BasisSetSize = lcob->BasisSetSize;
bs_eta.setIdentity(false);
*(bs_eta.C) = *(lcob->C);

12
src/QMCWaveFunctions/tests/test_soa_cusp_corr.cpp
View File

@ -129,14 +129,14 @@ TEST_CASE("applyCuspInfo", "[wavefunction]")
REQUIRE(lcob != NULL);
LCAOrbitalSet phi = LCAOrbitalSet(lcob->myBasisSet);
phi.setOrbitalSetSize(lcob->OrbitalSetSize);
phi.BasisSetSize = lcob->BasisSetSize;
LCAOrbitalSet phi = LCAOrbitalSet(lcob->myBasisSet, lcob->isOptimizable());
phi.setOrbitalSetSize(lcob->getOrbitalSetSize());
phi.BasisSetSize = lcob->getBasisSetSize();
phi.setIdentity(false);
LCAOrbitalSet eta = LCAOrbitalSet(lcob->myBasisSet);
eta.setOrbitalSetSize(lcob->OrbitalSetSize);
eta.BasisSetSize = lcob->BasisSetSize;
LCAOrbitalSet eta = LCAOrbitalSet(lcob->myBasisSet, lcob->isOptimizable());
eta.setOrbitalSetSize(lcob->getOrbitalSetSize());
eta.BasisSetSize = lcob->getBasisSetSize();
eta.setIdentity(false);
*(eta.C) = *(lcob->C);

6
tests/molecules/He_ae/CMakeLists.txt
View File

@ -1,6 +1,8 @@
# Deterministic test for QMC optimizer
IF (NOT MIXED_PRECISION)
IF (MIXED_PRECISION OR (QMC_COMPLEX AND NOT ENABLE_SOA))
MESSAGE_VERBOSE("Skipping optimizer deterministic He tests in mixed precision (QMC_MIXED_PRECISION=1) or complex with aos (QMC_COMPLEX=1 and ENABLE_SOA=0)")
ELSE()
SIMPLE_RUN_AND_CHECK(
deterministic-He_ae-opt
"${CMAKE_SOURCE_DIR}/tests/molecules/He_ae"
@ -19,7 +21,5 @@ IF (NOT MIXED_PRECISION)
TRUE
2 He_VMC_SCALARS # VMC
)
ELSE()
MESSAGE_VERBOSE("Skipping optimizer deterministic tests in mixed precision (QMC_MIXED_PRECISION=1)")
ENDIF()

16
tests/test_automation/jenkins_oxygen_cpu.sh
View File

@ -29,10 +29,12 @@ mkdir -p /dev/shm/${BUILD_TAG}-build
cd /dev/shm/${BUILD_TAG}-build
cmake -DQMC_COMPLEX=0 -DQMC_MIXED_PRECISION=0 -DENABLE_SOA=0 -DBUILD_AFQMC=1 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}
@ -51,10 +53,12 @@ mkdir -p ${BUILD_TAG}-build
cd ${BUILD_TAG}-build
cmake -DQMC_COMPLEX=0 -DQMC_MIXED_PRECISION=1 -DBUILD_AFQMC=1 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}
@ -73,10 +77,12 @@ mkdir -p ${BUILD_TAG}-build
cd ${BUILD_TAG}-build
cmake -DQMC_COMPLEX=1 -DQMC_MIXED_PRECISION=0 -DENABLE_SOA=0 -DBUILD_AFQMC=1 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}
@ -95,10 +101,12 @@ mkdir -p ${BUILD_TAG}-build
cd ${BUILD_TAG}-build
cmake -DQMC_COMPLEX=1 -DQMC_MIXED_PRECISION=1 -DBUILD_AFQMC=1 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}

16
tests/test_automation/jenkins_oxygen_cuda.sh
View File

@ -30,10 +30,12 @@ mkdir -p /dev/shm/${BUILD_TAG}-build
cd /dev/shm/${BUILD_TAG}-build
time cmake -DQMC_COMPLEX=0 -DQMC_MIXED_PRECISION=0 -DENABLE_SOA=0 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" -DENABLE_CUDA=1 /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}
@ -52,10 +54,12 @@ mkdir -p ${BUILD_TAG}-build
cd ${BUILD_TAG}-build
time cmake -DQMC_COMPLEX=0 -DQMC_MIXED_PRECISION=1 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" -DENABLE_CUDA=1 /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}
@ -74,10 +78,12 @@ mkdir -p ${BUILD_TAG}-build
cd ${BUILD_TAG}-build
time cmake -DQMC_COMPLEX=1 -DQMC_MIXED_PRECISION=0 -DENABLE_SOA=0 -DBUILD_AFQMC=1 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" -DENABLE_CUDA=1 /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}
@ -96,10 +102,12 @@ mkdir -p ${BUILD_TAG}-build
cd ${BUILD_TAG}-build
time cmake -DQMC_COMPLEX=1 -DQMC_MIXED_PRECISION=1 -DBUILD_AFQMC=1 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" -DENABLE_CUDA=1 /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}

16
tests/test_automation/jenkins_oxygen_gpu.sh
View File

@ -30,10 +30,12 @@ mkdir -p /dev/shm/${BUILD_TAG}-build
cd /dev/shm/${BUILD_TAG}-build
time cmake -DQMC_COMPLEX=0 -DQMC_MIXED_PRECISION=0 -DENABLE_SOA=0 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" -DQMC_CUDA=1 /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}
@ -52,10 +54,12 @@ mkdir -p ${BUILD_TAG}-build
cd ${BUILD_TAG}-build
time cmake -DQMC_COMPLEX=0 -DQMC_MIXED_PRECISION=1 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" -DQMC_CUDA=1 /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}
@ -74,10 +78,12 @@ mkdir -p ${BUILD_TAG}-build
cd ${BUILD_TAG}-build
time cmake -DQMC_COMPLEX=1 -DQMC_MIXED_PRECISION=0 -DENABLE_SOA=0 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" -DQMC_CUDA=1 /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}
@ -96,10 +102,12 @@ mkdir -p ${BUILD_TAG}-build
cd ${BUILD_TAG}-build
time cmake -DQMC_COMPLEX=1 -DQMC_MIXED_PRECISION=1 -DCMAKE_C_COMPILER="mpicc" -DCMAKE_CXX_COMPILER="mpicxx" -DQMC_CUDA=1 /dev/shm/${BUILD_TAG}-src 2>&1 | tee cmake.out
if [[ $? -ne 0 ]] ; then exit 1 ; fi
make -j 8
ctest -L unit --output-on-failure
if [[ $? -ne 0 ]] ; then exit 1 ; fi
ctest -L unit --output-on-failure
ret=$?
if [[ ${ret} -ne 0 ]] ; then
exit_code=${ret}