indexedTensor_tensor_solve.cpp 2.66 KB
Newer Older
Baum's avatar
Baum committed
1
// Xerus - A General Purpose Tensor Library
Sebastian Wolf's avatar
Sebastian Wolf committed
2
// Copyright (C) 2014-2016 Benjamin Huber and Sebastian Wolf. 
Baum's avatar
Baum committed
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
// 
// Xerus is free software: you can redistribute it and/or modify
// it under the terms of the GNU Affero General Public License as published
// by the Free Software Foundation, either version 3 of the License,
// or (at your option) any later version.
// 
// Xerus is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU Affero General Public License for more details.
// 
// You should have received a copy of the GNU Affero General Public License
// along with Xerus. If not, see <http://www.gnu.org/licenses/>.
//
// For further information on Xerus visit https://libXerus.org 
// or contact us at contact@libXerus.org.

20
/**
21
22
23
* @file
* @brief Implementation of the indexed tensor / operator.
*/
24

Sebastian Wolf's avatar
Sebastian Wolf committed
25
#include <xerus/misc/internal.h>
Sebastian Wolf's avatar
Sebastian Wolf committed
26
27
#include <xerus/misc/containerSupport.h>

28
#include <xerus/index.h>
29
30
#include <xerus/indexedTensor.h>
#include <xerus/indexedTensorMoveable.h>
Sebastian Wolf's avatar
Sebastian Wolf committed
31
#include <xerus/tensor.h>
32

Baum's avatar
Baum committed
33
namespace xerus {
34
35
	internal::IndexedTensorMoveable<Tensor> operator/ (internal::IndexedTensorReadOnly<Tensor>&& _b, internal::IndexedTensorReadOnly<Tensor>&& _A) {
		_A.assign_indices();
36
		_b.assign_indices();
Sebastian Wolf's avatar
Sebastian Wolf committed
37
		
Sebastian Wolf's avatar
Sebastian Wolf committed
38
		size_t extraDims = 0;
39
40
41
42
		std::vector<Index> orderA;
		std::vector<Index> orderB;
		std::vector<Index> orderX;
		
Sebastian Wolf's avatar
Sebastian Wolf committed
43
		// If possible we don't want to reorder A, so first divide A into those shared with b and x. 
44
		for(const Index& idx : _A.indices) {
45
46
47
48
49
50
51
			if(misc::contains(_b.indices, idx)) {
				orderA.push_back(idx);
			} else {
				orderX.push_back(idx);
			}
		}
		
Sebastian Wolf's avatar
Sebastian Wolf committed
52
53
		// This far the dimensions of A and B coincide
		orderB = orderA;
54
		
Sebastian Wolf's avatar
Sebastian Wolf committed
55
56
		// Add the second part of indices in the order obtained for X
		orderA.insert(orderA.end(), orderX.begin(), orderX.end());
57
		
Sebastian Wolf's avatar
Sebastian Wolf committed
58
59
		// Now complete indices of b and x with those not shared with A ( in order of b as we don't want to reorder b if possible).
		for(const Index& idx : _b.indices) {
60
			if(!misc::contains(_A.indices, idx)) {
Sebastian Wolf's avatar
Sebastian Wolf committed
61
62
63
64
65
				orderB.push_back(idx);
				orderX.push_back(idx);
				for(size_t i = 0; i < idx.span; ++i) {
					extraDims++;
				}
66
			}
67
68
		}
		
69
		// If indices coincide no reordering occours (only shared data pointer is set).
Sebastian Wolf's avatar
Sebastian Wolf committed
70
		Tensor reorderedA, reorderedB;
71
		reorderedA(orderA) = std::move(_A);
Sebastian Wolf's avatar
Sebastian Wolf committed
72
		reorderedB(orderB) = std::move(_b);
73
		
74
		internal::IndexedTensorMoveable<Tensor> tmpX(new Tensor(), std::move(orderX));
Sebastian Wolf's avatar
Sebastian Wolf committed
75
		
76
77
		//solve_least_squares(*tmpX.tensorObject, reorderedA, reorderedB, extraDims);
		solve(*tmpX.tensorObject, reorderedA, reorderedB, extraDims);
78
79
80
		
		return tmpX;
	}
81
} // namespace xerus