diff --git a/tests/test_stream.py b/tests/test_stream.py
index c373e26a..66d63ad3 100644
--- a/tests/test_stream.py
+++ b/tests/test_stream.py
@@ -373,6 +373,20 @@ def test_mixing_balance():
     ms['l'].mix_from(streams)
     assert ms['l'].imol['Water'] == total_flow
 
+def test_mixing_phases():
+    tmo.settings.set_thermo(['Water'], cache=True)
+    ms = tmo.MultiStream(None, l=[('Water', 1)], g=[('Water', 2)])
+    stream = tmo.Stream(Water=1, phase='L')
+    ms.mix_from([ms, stream], conserve_phases=True)
+    assert ms.phases == ('L', 'g', 'l') # New phase is created
+    assert ms.imol['L', 'Water'] == ms.imol['l', 'Water'] == 1
+    
+    ms = tmo.MultiStream(None, l=[('Water', 1)], g=[('Water', 2)])
+    stream = tmo.Stream(Water=1, phase='L')
+    ms.mix_from([ms, stream]) 
+    assert ms.phases == ('g', 'l') # 'L' phase gets added to 'l' phase
+    assert ms.imol['l', 'Water'] == 2
+
 def test_mixing_pressure():
     P = 111458
     s1 = tmo.Stream(None, Water=1000, P=P)
@@ -431,4 +445,5 @@ def test_vle_critical_pure_component():
     test_vle_critical_pure_component()
     test_critical()
     test_mixture()
+    test_mixing_phases()
     test_mixing_pressure()
\ No newline at end of file
diff --git a/thermosteam/_phase.py b/thermosteam/_phase.py
index eed9de35..ffff55e5 100644
--- a/thermosteam/_phase.py
+++ b/thermosteam/_phase.py
@@ -92,7 +92,7 @@ class PhaseIndexer:
     True
     
     """
-    __slots__ = ('_index',)
+    __slots__ = ('_index', '_phases', '_compatibility')
     _index_cache = {}
     
     def __new__(cls, phases):
@@ -102,32 +102,40 @@ def __new__(cls, phases):
             self = cache[phases]
         else:
             cache[phases] = self = new(cls)
-            self._index = index = {j:i for i,j in enumerate(sorted(phases))}
-            index[...] = slice(None) 
+            self._phases = phases = tuple(sorted(phases))
+            self._index = index = {j: i for i,j in enumerate(phases)}
+            for phase, n in tuple(index.items()):
+                if phase.isupper():
+                    phase = phase.lower()
+                else:
+                    phase = phase.upper()
+                if phase not in index: index[phase] = n
+            self._compatibility = ''.join([i.lower() for i in phases])
+            index[...] = slice(None)
         return self
     
+    @property
+    def phases(self):
+        return self._phases
+    
+    def __contains__(self, phase):
+        return phase in self._index
+    
+    def compatible_with(self, other):
+        return self._compatibility == other._compatibility
+    
     def __call__(self, phase):
         try:
             return self._index[phase]
         except:
-            if phase.isupper():
-                phase = phase.lower()
-            else:
-                phase = phase.upper()
-            try:
-                return self._index[phase]
-            except:
-                raise UndefinedPhase(phase)        
-    
-    @property
-    def phases(self):
-        return tuple(self._index)[:-1]
+            raise UndefinedPhase(phase)        
     
     def __reduce__(self):
-        return PhaseIndexer, (self.phases,)
+        return PhaseIndexer, (self._phases,)
     
     def __repr__(self):
-        return f"{type(self).__name__}({list(self.phases)})"
+        return f"{type(self).__name__}({[*self._phases]})"
+
 
 class Phase:
     __slots__ = ('_phase',)
diff --git a/thermosteam/_stream.py b/thermosteam/_stream.py
index c1701f6a..e991cfcd 100644
--- a/thermosteam/_stream.py
+++ b/thermosteam/_stream.py
@@ -1445,9 +1445,9 @@ def mix_from(self, others, energy_balance=True, vle=False, Q=0., conserve_phases
                     self.vle(T=self.T, P=P)
                 self.reduce_phases()
             else:
-                if energy_balance: H = sum([i.H for i in streams], Q)
-                self._imol.mix_from([i._imol for i in streams])
-                if energy_balance and not self.isempty():
+                if energy_balance: 
+                    self._imol.mix_from([i._imol for i in streams])
+                    H = sum([i.H for i in streams], Q)
                     if conserve_phases: 
                         self.H = H
                     else:
@@ -1457,6 +1457,8 @@ def mix_from(self, others, energy_balance=True, vle=False, Q=0., conserve_phases
                             self.phases = self.phase + ''.join([i.phase for i in others])
                             self._imol.mix_from([i._imol for i in streams])
                             self.H = H
+                else:
+                    self._imol.mix_from([i._imol for i in streams])
                 
     def split_to(self, s1, s2, split, energy_balance=True):
         """
diff --git a/thermosteam/indexer.py b/thermosteam/indexer.py
index 4665a5c7..d6deb2ea 100644
--- a/thermosteam/indexer.py
+++ b/thermosteam/indexer.py
@@ -774,11 +774,13 @@ def sum_across_phases(self):
     
     def copy_like(self, other):
         if self is other: return
+        phase_indexer = self._phase_indexer
         if isinstance(other, ChemicalIndexer):
             self.empty()
             other_data = other.data
-            self._expand_phases(other.phase)
-            phase_index = self.get_phase_index(other.phase)
+            phase = other.phase
+            if phase not in phase_indexer: self._expand_phases(phase)
+            phase_index = phase_indexer(phase)
             if self.chemicals is other.chemicals:
                 self.data.rows[phase_index].copy_like(other_data)
             else:
@@ -786,17 +788,32 @@ def copy_like(self, other):
                 left_index, right_index = index_overlap(self._chemicals, other._chemicals, [*other_data.nonzero_keys()])
                 self.data.rows[phase_index][left_index] = other_data[right_index] 
         else:
+            other_phase_indexer = other._phase_indexer
             if self.chemicals is other.chemicals:
-                self._expand_phases(other.phases)
-                self.data.copy_like(other.data)
+                if phase_indexer is other_phase_indexer:
+                    self.data.copy_like(other.data)
+                elif phase_indexer.compatible_with(other_phase_indexer):
+                    self.empty()
+                    data = self.data
+                    for i, j in other: data[phase_indexer(i)] = j
+                else:
+                    self._expand_phases(other._phases)
+                    self.data.copy_like(other.data)
             else:
                 self.empty()
-                self._expand_phases(other.phases)
                 other_data = other.data
+                data = self.data
                 left_index, other_data = index_overlap(self._chemicals, other._chemicals, [*other_data.nonzero_keys()])
-                self.data[:, left_index] = other_data[:, right_index]
+                if phase_indexer is other_phase_indexer:
+                    data[:, left_index] = other_data[:, right_index]
+                elif phase_indexer.compatible_with(other_phase_indexer):
+                    for i, j in other: data[phase_indexer(i)] += j
+                else:
+                    self._expand_phases(other._phases)
+                    data[:, left_index] = other_data[:, right_index]
+                    
     
-    def _expand_phases(self, other_phases):
+    def _expand_phases(self, other_phases=None):
         phases = self._phases
         other_phases = set(other_phases)
         new_phases = other_phases.difference(phases)
@@ -822,10 +839,22 @@ def mix_from(self, others):
             else: raise ValueError("can only mix from chemical or material indexers")
         other_phases = [i.phase for i in chemical_indexers]
         for i in material_indexers: other_phases.extend(i._phases)
-        self._expand_phases(other_phases)
+        other_phases = set(other_phases)
+        phase_indexer = self._phase_indexer
+        new_phases = [i for i in other_phases if i not in phase_indexer]
         phases = self._phases
+        if new_phases: self._expand_phases(other_phases)
         scp_data = {i: [] for i in phases} # Same chemicals by phase
         dcp_data = {i: [] for i in phases} # Different chemicals by phase
+        for i in other_phases.difference(phases):
+            if i.isupper():
+                ilow = i.lower()
+                scp_data[i] = scp_data[ilow]
+                dcp_data[i] = dcp_data[ilow]
+            else:
+                iup = i.upper()
+                scp_data[i] = scp_data[iup]
+                dcp_data[i] = dcp_data[iup]
         for i in material_indexers:
             ichemicals = i._chemicals
             idata = i.data
@@ -835,7 +864,7 @@ def mix_from(self, others):
             else:
                 left_index, right_index = index_overlap(chemicals, ichemicals, idata.nonzero_keys())
                 for i, j in zip(i._phases, i.data.rows):
-                    dcp_data[i].append((idata, left_index, right_index))
+                    dcp_data[i].append((j, left_index, right_index))
         for i in chemical_indexers:
             ichemicals = i._chemicals
             idata = i.data